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Abstract. In the series of 3d (tig) 1 perovskites, SrV03-CaV0 3 -LaTi0 3 -YTi03, 
the transition-metal d electron becomes increasingly localized and undergoes a 
Mott transition between CaVC>3 and LaTi03. By defining a low-energy Hubbard 
Hamiltonian in the basis of Wannier functions for the ti g LDA band and solving it 
in the single-site dynamical mean-field (DMFT) approximation, it was recently shown 
(Pavarini et al (2004) Phys. Rev. Lett. 92 176403) that simultaneously with the Mott 
transition there occurs a strong suppression of orbital fluctuations due to splitting of 
the tig levels. The present paper reviews and expands this work, in particular in the 
direction of exposing the underlying chemical mechanisms by means of ab initio LDA 
Wannier functions generated with the Ath order muffin-tin orbital (NMTO) method. 
The Wannier functions for the occupied oxygen-p band illustrate the importance of 
oxygen-p to large cation-d covalency for the progressive GdFeC>3-type distortion along 
the series. The oxygen-p orbitals which pder-bond to the cations are the same as those 
which pdir-bond to the transition-metal tig orbitals. As a consequence, the Wannier 
functions for the ti g band exhibit residual covalency between the transition- metal t% g , 
the large cation-d, and the oxygen-p states. This residual covalency, which increases 
along the series, turns out to be responsible not only for the splittings, A, of the ti g 
levels, but also for non-cubic perturbations of the hopping integrals, both of which 
are decisive for the Mott transition. We find good agreement with the optical and 
photoemission spectra for all four materials, with the crystal-field splittings and orbital 
polarizations recently measured for the titanates, and with the metallization volume 
(pressure) for LaTiC>3. The metallization volume for YTiC>3 is predicted and the role of 
the Jahn- Teller distortion is discussed. For use in future many-body calculations, we 
tabulate the ti g on-site and hopping matrix elements for all four materials and give an 
analytical expression for the orthorhombic Hamiltonian in the k + Q representation. 
Using conventional super-exchange theory, our on-site and hopping matrix elements 
reproduce the observed magnetic orders in LaTiC>3 and YTiC>3, but the results are 
sensitive to detail, in particular for YTi03 where, without the JT distortion, the 
magnetic order would be antiferromagnetic C- or A-type, rather than ferromagnetic. 
It is decisive that upon increasing the GdFe03-type distortion, the nearest-neighbour 
hopping between the lowest and the upper-level Wannier functions becomes stronger 
than the hopping between the lowest-level Wannier functions. Finally, we show that the 
non-cubic perturbations responsible for this behaviour make it possible to unfold the 
orthorhombic ti g LDA bandstructure to a pseudo-cubic zone. In this zone, the lowest 
band is separated from the two others by a direct gap and has a width, Wj, which is 
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significantly smaller than that, W, of the entire t^g band. The progressive GdFeC>3- 
type distortion thus favours electron localization by decreasing W, by increasing A/W, 
and by decreasing Wi/W. Our conclusions concerning the roles of GdFeC>3-type and 
JT distortions agree with those of Mochizuki and Imada (2003 Phys. Rev. Lett. 91 
167203). 
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1. Introduction 

Transition-metal perovskites have been studied for half a century, and most intensively 
during the last decade, for their fascinating electronic and magnetic properties arising 
from narrow 3d bands and strong Coulomb correlations [1, 2, 3]. The 3d(t 2g ) 1 series 
SrV03-CaV0 3 -LaTi0 3 -YTi03 is a paradigm because it has no complicating multiplet 
effects, a progressing structural distortion illustrated at the top of figure 1, and greatly 
varying electronic properties: while Sr and Ca vanadate are correlated metals, with 
optical mass enhancements of respectively ~3 and ~4 [4, 5, 6], La and Y titanate are 
Mott insulators, with gaps of respectively 0.2 and 1 eV [7]. These two Mott insulators, 
which essentially mark the end points of a series of rare-earth titanates RTi03,[3, 8] have 
very different metallization pressures, 11 GPa for LaTi03 and much larger for YTiC>3 [9]. 
Moreover, they exhibit different orbital physics [10] and, at low temperature, LaTi0 3 is 
a 3-dimensional (G-type) antiferromagnet with T/v=150K [11] and a small moment of 
0.57 he [12], while YTi0 3 is a ferromagnet with a low Curie temperature of T c =S0K 
and a good-sized moment of 0.8/x B [11, 13]. 

In the Mott-Hubbard picture the metal-insulator transition occurs when the ratio 
of the on-site Coulomb repulsion to the one-electron bandwidth exceeds a critical value, 
(U/W) c . As figure 1 shows, in the ABO3 perovskites the B 3d ions are on a nearly 
cubic (orthorhombic) lattice and at the centres of corner-sharing 06 octahedra. The 
3d band thus splits into covalent Op Bd n antibonding t 2g bands and covalent O p 
B d a antibonding e g bands, of which the former lie lower, have less O character, and 
couple less to the octahedra than the latter. Simple theories for the d 1 perovskites [2] 
are based on a Hubbard model with three independent, two-dimensional, degenerate, ~ 
filled tig bands per B ion, and the variation of the electronic properties along the series 
is ascribed to a progressive reduction of W due to the increased bending of the pdir 
hopping paths, the B-O-B bonds seen in figure 1. 

But this need not be the whole story, because the value of (U/W) c is expected 
[14] to decrease with decreasing degeneracy. Such a decrease of degeneracy can be 
achieved by splitting the t 2g levels by merely ZW, the reduced bandwidth associated 
with quasiparticle excitations in the correlated metal [15]. As a consequence, if the level- 
splitting increased along the series, this could significantly influence the Mott transition. 

However, unlike in e g band perovskites, where large (10%) cooperative Jahn- Teller 
(JT) distortions of the oxygen octahedron indicate that the orbitals are not degenerate 
but spatially ordered, the octahedron in the t 2g band perovskites is nearly perfect. For 
that reason the t 2g orbitals have often been assumed to be degenerate. If that is the 
case, quantum fluctuations will lead to an orbital liquid [10, 16, 17] rather than orbital 
ordering in the Mott insulating phase. The observation of an isotropic, small gap spin- 
wave spectrum both in antiferromagnetic LaTi0 3 and also in ferromagnetic YTi0 3 [10] 
has lent support to this orbital-liquid scenario, because if the orbital moments were 
quenched, such a spectrum would seem accidental. On the other hand, the predicted 
contribution to the specific heat from the orbital liquid has not been observed in LaTi0 3 
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Figure 1. Crystal structures and electronic bandstructures for the series of 3d(t2 g ) 
orthorhombic ABO3 perovskites considered in this paper; A (green) B (red) O (blue). 
The bottom part shows densities of one-electron states (DOSs) calculated in the LDA 
for the real structures (right-hand panels) and for hypothetical, cubic structures with 
the same volumes (left-hand panels). The green, red, and blue DOSs are projected 
onto, respectively, Ad, B3d, and 2p orthonormal orbitals [41]. The B3d(i2 9 ) bands 
are positioned around the Fermi level (zero of energy) and their widths, W, decrease 
from ~3 to ^2eV along the series. The much wider B3d(e g ) bands are at higher 
energies. This figure resulted from linear muffin-tin orbitals (LMTO) calculations in 
which the energies, e v m, of the linear, partial- wave expansions were chosen at the 
centres of gravity of the occupied, partial DOS. Since those energies are in the O p 
band, the LMTO errors proportional to (e — f. v m) slightly distort the -unoccupied 
parts of the DOS. 

[18]. Moreover, a 3% JT stretch of one of the basal O squares into a rectangle was 
recently discovered in LaTi03 [12, 19]. This is of similar magnitude as the JT distortion 
known to exist in YTiC>3 where, however, the square is stretched into a rhomb [20]. 

By ab initio calculation of the Wannier functions of the LDA t 2g band, it was 
recently found [21] that in LaTiOs and YTiC>3 the t<i g degeneracy is lifted at the classical 
level. This is not due to the small JT distortions, but to the GdFe0 3 -type distortion 
which tilts and rotates the corner-sharing octahedra as illustrated in figure 2. As we 
shall see in the present paper, this distortion is partly driven by the covalency between 
occupied oxygen p states and empty A-cation d states, which pulls each 01 (02) closer 
to one (two) of its four nearest A neighbours. As a result, each A cation has 4 of its 
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12 near oxygens pulled closer in. In addition, the A cube gets distorted so that one 
diagonal becomes the shortest. The t 2g degeneracy is now lifted, essentially by residual 
covalent interactions between empty Ad orbitals and full Bt 2g and O p orbitals, the 
details of which are different in LaTiOs and YTi03. These residual covalent interactions 
not only perturb the on-site, but also the hopping matrix elements of the Hamiltonian. 

Already forty years ago, Goodenough [22] speculated that covalency between 
occupied Op and empty Ad orbitals may be a driving force behind the GdFe0 3 - 
type distortion, and this hypothesis was recently supported in an extensive series of 
semiempirical simulations [23]. 

That the GdFeO"3-type distortion lifts the t 2g degeneracy in the titanates had been 
realized slightly earlier by Mochizuki and Imada [24], but only after a long search by 
them and other groups for a model Hamiltonian which could reproduce all observed 
magnetic and orbital orderings [25, 26, 27, 28, 29]. These model Hartree-Fock and 
strong-coupling studies for the entire family of 3d n perovskites mapped out the roles 
played JT distortions, spin-orbit coupling, e g degrees of freedom, orbital misalignment 
caused by the GdFe0 3 distortion, and, finally, by the electrostatic field and hybridization 
from the A cations which enter via the GdFe0 3 -type distortion. 

What enabled the single, independent study in reference [21] to reach the 
same conclusion concerning the role of the GdFe03-type distortion, was the use of 
parameter-free density-functional (LDA) theory to generate t 2g Wannier functions for a 
representative series of real materials by means of a new technique [30]. The point is that 
it is virtually impossible to know a priori which on-site and hopping matrix elements are 
important, and what their values are. To extract them with the required accuracy from 
experiments and/or from LDA bandstructures is often impossible, and in this respect 
the t 2g perovskites are particularly nasty: at first sight, the t 2g Hamiltonian of these 
nearly cubic materials is the simplest possible, but cation covalency in the presence of 
GdFe03-type distortion makes it far more complicated than, e.g., the e g Hamiltonian 
relevant for the t\ g e l g manganites, which have similar GdFe03 distortions. The reason 
is that the t 2g orbitals interact with the same oxygen p orbitals as the cation orbitals 
do. This is not the case for e g orbitals. 

Now the problem that relevant parameters might be overlooked in model 
calculations occurs less frequently with the LDA+?7 method, because this method uses a 
complete basis set. In fact, at about the time when Mizokawa and Fujimori [25] carried 
out their pioneering model calculations, Solovyev, Hamada, Sawada, and Terakura 
[31, 32] performed LDA+£7 calculations in which the value of the on-site Coulomb 
repulsion, U, was adjusted to the observed optical gaps. Those calculations yielded 
the correct magnetic orders in LaTi03 and YTi03, but the underlying mechanism was 
not recognized. The magnetic moment and orbital order predicted for YTi03 were 
subsequently confirmed by NMR [33] and neutron scattering [34]. 

The static mean-field approximation used in Hartree-Fock and LDA+[7 calculations 
cannot describe the paramagnetic-metal to paramagnetic-insulator (Mott) transition. 
To be more specific, this approximation is inadequate for strongly correlated metals 
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Figure 2. Crystal structure of a GdFe03-distorted ABO3 perovskite (YTiOa) in Pbnm 
symmetry. The horizontal AO-planes perpendicular to the z-axis are mirrors. The 
global x-axis points outward in the front and the y-axis to the right. The B-sites (red) 
form a primitive monoclinic, nearly cubic lattice with the following translation vectors: 
R* = [1 + [ft + a) /2] x + [{ft - a) /2] y, R y = [(ft - a) /2] x + [1 + (ft + a) /2] y, and 
R 2 = (1 + 7) z with a, ft and 7 small. The latter vector is orthogonal to the two former, 
which have the same length, but need not be quite orthogonal to each other. For 
LaTi0 3 (YT1O3), |x| = \y\ = |z| = 397.1 (385.8) pm and a = 33 (-258), = 1 (407), 
7 = —33 (—137) x 10~ 4 . The position of a B-site is: R = xR x + yH y + zR z with 
x, y, and z integers. The structure is orthorhombic with 4 ABO3 units per cell, e.g., 
the four in the front plane, which we shall label 1 (bottom left, xyz=000), 2 (bottom 
right, 010), 3 (top left, 001) and 4 (top right, 011). The orthorhombic translation 
vectors are: a = (x — y) (1 + a) — R x — Ry, b = (x + y) (1 + ft) = R x + R y , and 
c = 2z (1 + 7) = 2H Z . Reflection in a vertical 6c-plane containing a B-site (a ^ —a, or 
equivalently, x <-> y), followed by a translation R. y = | (b + a) , takes the crystal into 
itself (glide plane). All A-ions (yellow) and all B-ions (red) are equivalent, but there 
are two kinds of oxygen (blue): 01 in a mirror A-plane, and 02 in a vertical, buckled 
A-plane. Proceeding along the series in figure 1, the GdFe03-type distortion tilts the 
corner-sharing octahedra by 0, 9, 13 (12), and 20° around the 6-axes in alternating 
directions, and rotates them around the c-axis by 0, 7, 10 (9), and 13° in alternating 
directions [12, 20, 42, 43]. Here, the values in parentheses are from the older data [44]. 
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and for Mott insulators at temperatures above the magnetic ordering temperature. For 
such purposes, the dynamical mean-field approximation (DMFT) [35] has recently been 
developed and applied to solve low-energy Hubbard Hamiltonians derived from the LDA 
[36]. The calculations reported in reference [21] employed the above-mentioned basis of 
LDA t2 g Wannier functions, and a single, adjusted value of U for all four materials. This 
implementation of the LDA+DMFT approach properly describes orbital fluctuations 
in the multiband Hubbard model by including the off -diagonal mm'-matrix elements 
of the self-energy matrix, taken to be local in DMFT. All earlier implementations 
had used a scalar self-energy, as is appropriate for cubic systems [37, 38]. A recent 
LDA+DMFT calculation for Lai-^Sr^TiOa used the eigenrepresentation of the on-site 
LDA Hamiltonian and then neglected the off-diagonal elements of the self-energy [39], 
apparently a reasonable approximation for this system. 

This paper is a pedagogical review of the calculations reported in reference [21] and 
a presentation of many new results. The use of LDA Wannier functions to bring out the 
materials aspects is emphasized. The paper takes off (section 2) with a discussion 
of the chemistry: after analyzing the high-energy part of the LDA bandstructures 
displayed in figure 1 , we demonstrate the role of O-A covalency in driving the progressive 
GdFe0 3 distortion. In particular, we visualize the bonds by means of the O 2p Wannier 
functions. In section 3, we zoom in on the low-energy LDA t 2 g-bands and discuss 
their Wannier functions, first in the cubic xy, yz, xz-representation and then in the 
crystal-field representation. The influence of the GdFeO"3-type and JT distortions on 
the orbital energies, inter-orbital couplings, and the bandstructures is discussed in 
detail, and the on-site and hopping integrals are tabulated for use in future many- 
body calculations. In section 4 we set up the t2 9 Hubbard Hamiltonian and explain how 
it is solved with the DMFT many-body technique. Using the same value of U for all four 
materials, the resulting high-temperature electronic structures are presented in section 
5. They reproduce the increased localization observed along the series, including the 
mass enhancements, the Mott transition, and the gap sizes. For the Mott insulators, 
it turns out that the Coulomb correlation localizes the electron almost exclusively in 
one orbital, and more so for YTiO"3 than for LaTiO"3, and that this orbital is the 
eigenfunction of the lowest t2 9 level in the LDA. The roles of band-width reduction 
and level splittings are analyzed, and it is concluded that these alone do not suffice to 
explain the strong decrease of (U/W) c found by DMFT when progressing through the 
series. We then compute the metallization volumes for the titanates, compare with most 
recent high-pressure experiments [9], and discuss the role of the GdFe0 3 -type and JT 
distortions. In order to further elucidate the difference in the electronic structure of the 
Mott insulators LaTiOa and YTi0 3 we present calculations of the onset of the optical 
conductivity. In section 6 we use our computed orbital orders and hopping integrals 
to calculate the magnetic exchange couplings, J se , within conventional super-exchange 
theory. Reasonable agreement with the observed magnetic orders is obtained, but -as 
evidenced by the extensive model calculations [25, 26, 27, 28, 29], and as pointed out 
by Ulrich et al. [40] - such J se values are extremely sensitive to detail, in particular in 
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Figure 3. LDA bandstructure of cubic SrVC>3. The Brillouin zone (BZ) is shown in 
blue in figure 9. As is discussed in appendix A, the bands obtained with the downfolded, 
truly minimal V 3d fog) NMTO basis (red) are indistinguishable from those obtained 
with the full NMTO basis (black). 



YTiC>3 where we find the small JT distortion to be decisive for the ferromagnetism. In 
section 7 we uncover the ingredient so far missing in our understanding of the calculated 
trend in {U/W) c , namely the formation of a lowest subband, whose width is significantly 
smaller than W. Specifically, the iVth-order muffin-tin orbitals (NMTO) method enables 
us to show that the complicated orthorhombic LDA t2 9 bandstructures for the titanates 
can be approximately folded out to a pseudo-cubic zone and that, in this representation, 
the lowest band is separated from the two upper bands by a direct gap. The development 
of this orbitally ordered band from the three degenerate cubic bands is explained and 
the relation to the increasing tendency towards ferromagnetism with increasing GdFeC>3- 
distortion is pointed out. In section 8 we sum up our main conclusions. Appendix A 
explains how downfolding within the NMTO method [30] is used to construct truly 
minimal basis sets which pick out selected bands, such as the O p or B t 2g bands. When 
symmetrically orthonormalized, such a truly minimal basis set constitutes a set of atom- 
centred, highly localized Wannier functions. Finally, since a weak point of the present 
calculations is our use of the standard LMTO-ASA method [41] to generate the LDA 
potentials, we give the technical details in appendix B. 
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2. High-energy LDA bandstructures and O-A covalent mechanism of the 
GdFe0 3 -type distortion 

In the present section we shall use the LDA bandstructures shown in figures 1, 3, and 
4, as well as the oxygen p band Wannier functions, to demonstrate the role of O p — A 
d covalency for the GdFe0 3 -type distortion, and the role of O p — B d covalency for the 
stability of the octahedron and the splitting of the B 3d band into separate t 2g and e g 
bands. This will set the stage for understanding the role of covalency for the splitting 
of the t 2g levels. 

Apart from providing insights into the bonding, the LDA bands -except for the 
t 2g bands, which we shall treat separately with a Hubbard model- give information 
about the one-electron high-energy excitations. Specifically, we shall approximate these 
excitations by the LDA Hamiltonian plus a self-energy, E (e) , which only couples inside 
t 2g space and is independent of the crystal momentum, k. To compute this self-energy 
is the task of the LDA+DMFT to be considered in section 4. Future calculations 
might choose to include in the many-body calculation a larger basis set of Wannier 
functions than those describing merely the low-energy t 2g bands. For that reason, too, 
it makes sense to consider first the high-energy LDA band structure and the O p Wannier 
functions. 

2.1. Bandstructures 

The bottom part of figure 1 exhibits densities of states (DOSs) projected onto various 
groups of orbitals, Op, Bd, and Ad, of a standard set of nearly orthonormal LMTOs 
[41]. These DOS projections provide information about the mixing of characters due to 
hybridization between various kinds of orbitals. Consider for instance the panel relating 
to SrV0 3 : in the Op bands (mainly blue), we see V d character (red), and in the Yd 
bands (mainly red), we see O p character (blue). This hybridization between the Op 
and the Yd bands has pushed them apart, and since the oxygen bands are occupied 
and the V d bands are nearly empty, band-structure energy has been gained; this is O-B 
covalency. Some orbitals hybridize more than others, e.g., a bonds are stronger than 
7i bonds, and it is therefore the lower part of the Op band and the upper, e g part of 
the Vd band which have the most foreign character mixed in. The less familiar result 
of this figure is that it also exhibits a large amount of O-A (blue-green) covalency. We 
shall see that this is because each oxygen has two p orbitals a bonding with Sr, but only 
one a bonding with V. That, to some extent, compensates for the distance to Sr being 
a/2 longer than the distance to V. 

When we move along the series SrV03-CaV0 3 -LaTi0 3 -YTi03 a tilt is 
expected if the ionic radius is such that the Goldschmidt tolerance factor, 
(rA + ro) I [a/2 (t\b + ro)] , is smaller than 1. The radii satisfy: r Sr 2+ ~ r L a 3+ > r Ca 2+ ~ 
rya+, while rys+ < r Ti 4+ (in this scheme the d electron is taken as localized), and the 
tolerance factor decreases by about 10% along the series, although it is the same for 
CaV0 3 and LaTi0 3 . The progressive tilt is thus partly due to the fact that the size of 
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Figure 4. LDA bandstructure of orthorhombic LaTiC>3. The Brillouin zone is shown 
in red in figure 9. The bands obtained with the truly minimal (downfolded) O 2p 
NMTO basis (red) are indistinguishable from those obtained with the full NMTO basis 
(black). This is explained in appendix A. The basis functions of the truly minimal set 
are shown in figures 5 and 6. The recent structural data [12] was used. 

the A cation shrinks in relation to that of the BO3 octahedron. However, as pointed out 
by Woodward [23], this does not explain why the tilt is of the GdFeC>3-type. This type, 
he found by using empirical interatomic potentials and the extended Hiickel method, 
is unique in maximizing the O-A covalent bonding at the same time as minimizing the 
O-A repulsive overlap. Our LDA calculations support and detail the importance of O-A 
covalency for the GdFe0 3 -type distortion: 

It is well accepted that equilibrium crystal structures, such as those shown at the top 
of figure 1, may be computed ah initio with good accuracy using the LDA. In the bottom 
right-hand panels we therefore show the LDA densities of states (DOS) calculated for the 
real structures -determined experimentally though [20, 42, 43, 44]- and in the bottom 
left-hand panels we show the LDA DOS calculated for hypothetical, cubic structures 
with the same volume. Now, the energy gain associated with a structural distortion is 
approximately the gain in band-structure plus Madelung energy, so let us consider the 
trend in the former: 

For SrV03 the left- and right-hand panels are identical because the real structure 
of SrV03 is cubic. Each Sr ion is at the corner of a cube and has 12 nearest oxygens 
at the face centres. Going now to cubic CaV03 the empty 3d band of Ca lies lower 
and thereby closer to the oxygen 2p band than the empty 4d band of Sr. It is therefore 
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Figure 5. Basis functions of the truly minimal set of 2p NMTOs calculated for 
LaTiC>3 with the structure of reference [12]. Shown are the orbital shapes (constant- 
amplitude surfaces) with the ± signs labelled by red and blue. 01 is in a flat face of 
the distorted La cube and the 02s are in the buckled faces (see caption to figure 2). In 
column 1 we show the orbitals perpendicular to the faces (p z ) which exhibit symmetric 
Ti-O-Ti cr-bonds. In column 2 we show the orbitals in the top horizontal face (01 p x 
and p y ). Although they are equivalent to those in the bottom face shown in column 3, 
we show both for the sake of clarity. The 01 p x and p y orbitals, as well as those for 
02 shown in the following figure 6, exhibit weak, symmetric Ti-O-Ti 7r-bonds. Most 
importantly, however, the 01 p x and p y orbitals show asymmetric O-La cr-bonds. The 
latter, together with the 02 p x - La cr-bonds shown in figure 6, are responsible for the 
GdFe03-type distortion. This 0-A bonding is shown schematically in figure 7. 



conceivable that a GdFe03-type distortion which pulls some of the oxygen neighbours 
closer to the A ion and thereby increases the covalency with those, is energetically more 
favourable in CaVC>3 than in SrVC>3, and this is what the figure shows: an increase 
of the 2p-Ca3<i gap associated with the distortion in CaVC>3. The Ca3<i character is 
essentially swept out of the lower part of the V 3d band. 

When proceeding to the titanates, the A and B cations become 1st- rather than 3rd- 
nearest neighbours in the periodic table. The B 3c? band therefore moves up and the Ad 
band down with respect to the Op band. Hence, the 0-B covalency decreases and the O- 
A covalency increases. Most importantly, the A d band becomes nearly degenerate with 
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Figure 6. Same as figure 5, but for the oxygen orbitals (02 p x and p y ) in the buckled 
faces of the distorted La cube. 02 p x is seen to bind much stronger to La than 02 p y . 
The orbitals in the 1st and 3rd columns are equivalent, but for clarity we show both. 
Not shown are the orbitals on the rear face since they are equivalent to those in the 
front face shown in the 2nd column. The rear p x orbital binds strongly to La at the 
bottom right corner. The 0-A bonding is summarized in figure 7. 

the Ti3c? band, and more so for Y4d than for La 5c?. It is only the GdFe03-distortion 
which, through increase of the O 2p-A d hybridization, pushes the A d band above the 
Fermi level. This, as well as the concomitant lowering of the O 2p band with respect to 
€p, can be seen in the figure. 

In the series of trivalent rare-earth titanates, the GdFe03-type distortion -and the 
low-temperature magnetic properties- are known to change gradually from those of 
LaTiOa to essentially those of YTi03 [3]. This trend follows the decrease of the ionic 
radius and of the 5c? level position. 

2.2. Oxygen p bonds 

The covalent bonds can be visualized by a set of localized Wannier functions for the 
occupied bands, i.e. the O 2p bands. Here we shall not consider the much smaller 
contribution from the O 2s bands. In figure 4 we illustrate for the case of LaTi03 that 
we can construct a set of O p NMTOs which span the O 2p bands, and no other bands. 
This basis set contains as many orbitals as there are occupied bands and, hence, it is 
what we call a truly minimal basis set. Its inequivalent orbitals are shown in figures 5 
and 6. Each oxygen has two nearest B neighbours and four near A neighbours, and we 
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have chosen its three p orbitals such that one of them, call it p z , points approximately 
towards the B neighbours (see caption to figure 2), and the two other orbitals, labelled 
p x and p y , point approximately towards two of the four A neighbours. 

Only the central parts of these orbitals have respectively p x , p y , or p z character, 
because in order to describe the hybridization of the O 2p band with the B and A derived 
bands, the orbitals of a truly minimal set must have those characters folded into their 
tails. In fact, all partial- wave characters other than oxygen p are folded into the Op 
set. Now, as is explained in appendix A, the NMTOs are localized by construction, but 
they are not quite orthogonal. Therefore, in order to become a set of localized Wannier 
functions, the truly minimal set must be symmetrically orthonormalized. The sum of 
the squares of these Wannier functions is the valence charge density, and the sum of 
their orbital energies is the band-structure energy (except for a small contribution from 
the B 3d electron). Hence, the truly minimal set of O 2p NMTOs shown in figures 5 and 
6 visualizes the covalent bonds. 

From the first column in figure 5 we see that each oxygen forms one strong Ti-O-Ti 
covalent, symmetric a bond with its two Ti neighbours, the Ti character being 3d (e g ) 
with a bit of Asp. From the right-hand side of the figure, we see that each oxygen orbital 
in the horizontal, flat face of the distorted La cube (01 p x or p y ) forms a weak Ti-O-Ti 
covalent, symmetric 7r bond, the Ti character being predominantly 3d (£ 2ff ) • As may be 
seen from figure 6, the same holds for each oxygen orbital (02 p x or p y ) in the vertical, 
buckled faces. The respective O-Ti covalent interaction pushes the anti-bonding Ti 
3c? (e s )-like band well above the Ti 3d (i2 9 )-like band, as we saw in figure 1. 

Most importantly, however, figure 5 shows that each of the two oxygen orbitals in 
the flat faces (01 p x or p y ) forms an 0-A covalent, asymmetric a bond with one of the 
two La neighbours towards which it points. It is clearly seen how the weight of the 
orbital is shifted from one towards the other La ion. Moreover, we see that the bond 
with 01 p x is somewhat stronger than with 0\p y . The La character of the 01 p x orbital 
is 5g? 3x 2_ 1 with some by-mixing of 6sp. This d character is mostly t2 g , because with x 
and y along the cubic directions, it is d3^ x±y ~ ) 2_ 1 = ^^-d xy — ^d 3z 2_ 1 . Figure 6 shows 
that for the buckled faces of the distorted La cube, only one of the oxygen orbitals, 
02 p x , bonds significantly to La, and that this bonding is as strong as for 01 p x . As a 
result, we obtain the schematic picture of the 0-A covalent a bonds shown in figure 7. 
The resulting GdFe0 3 -type distortion shortens the two 01 bonds by respectively 17% 
(15%) and 11% (8%) of the average of the four Ol-La distances, and it shortens the 02 
bond by 16% (14%) of the average of the four 02-La distances [12]. Here the numbers 
in parentheses are from the older data [44]. For CaVOs the corresponding bond-length 
reductions are 10% and 4% for 01 and 12% for 02, while for YTi03, they are as large 
as 28% and 23% for 01, and 22% for 02. For YTi0 3 the shortest 0-Y distance is, in 
fact, only 10% longer than the O-Ti distance. 

Also indicated in figure 7 is the oxygen coordination of the A cations, which is 
reduced from 12-fold and cubic to 4-fold and roughly tetrahedral. The embedding of 
the unit in the entire structure can be understood by comparison with figure 2, where 
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Figure 7. Schematic representation of the O-A covalent bonds shown for LaTiC>3 in 
figures 5 and 6: 01 binds to two while 02 binds to one of the four A-neighbours. One 
of the two Ol-A bonds is relatively weak and is indicated by a short, red arrow. The 
resulting GdFe03-type distortion shortens the O-A bonds correspondingly. In CaV03, 
LaTi03, and YTi03, the shortest Ol-A bond is shortened by respectively 10, 17, and 
28% with respect to the average, the 2nd-shortest Ol-A bond by respectively 4, 11, 
and 23%, and the shortest 02-A bond by respectively 12, 16, and 22%. The oxygen 
coordination of the A-ion is reduced from 12 to 4, with two of the near oxygens being 
in the horizontal, flat face of the distorted A-cube, and the two others in one of the 
vertical, buckled faces. The A-B-A diagonal (orange bar) lying in the plane of the 
short, red arrows is shortened by respectively 3, 7, and 9% of the average. The unit 
shown is the front bottom left one (subcell 1) seen in e.g. figure 2. 



the unit is the one in the front bottom left corner. 
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Finally, the orange bar in figure 7 indicates that one of the four A-B-A diagonals, 
[111], is shortened by respectively 3, 7 (5), and 9% of the average A-B-A distance in 
CaVC>3, LaTiOs, and YTiC>3. The two A ions closest to B are those which bond weakest 
to the oxygen octahedron (red arrows in the figure). The corresponding distortion of 
the A cube is clearly visible in figure 2, and is presumably caused by the hard-core 
repulsion from the three nearest oxygens (the blue arrows in figure 7). Corresponding 
to this shortening of the [111] diagonal in subcell 1, is a stretching of the [111] diagonal 
by nearly the same amount. In YTi0 3 , the tilt of the oxygen octahedron around the 
b axis, i.e. towards the two Y atoms along [111], is particularly strong (20°, see figure 
2) and the corresponding A-B-A distance is shortened by 5%. None of the other three 
materials have such a second, short diagonal. In the distorted structure, each A ion 
(all equivalent) has two nearest B neighbours (two short diagonals). From the point of 
view of the A ion in the lower right corner of figure 7, i.e. the one at these B 

neighbours are at 110 and 111. 

In the following we shall turn to our primary interest, the development of the low- 
energy electronic structure along the series of d 1 perovskites. We shall see that the 
minute A character left behind in the B 3c? (t 2g ) band after most of this character has 
been swept away by O-A covalency and the concomitant GdFe0 3 -distortion, is decisive. 

3. t 2g Wannier functions and their Hamiltonian 

The physical properties of the d l perovskites are determined by the Wannier functions of 
the low-energy B 3d (t 2g ) bands, such as those shown in red in figure 3, and their on-site 
Coulomb repulsion. This, at least, is the working hypothesis of the present paper. 

In figure 8 we show the truly minimal t 2g NMTO basis set for the case of LaTi0 3 . 
In order to generate these orbitals, it is not necessary to choose local axes oriented after 
the oxygen octahedron: with x, y, and z referring to the global axes (see caption to 
figure 2) and the active channels specified simply as d xy , d yz , and d xz on each of the four 
B sites, each orbital automatically adjusts to its environment. In the present case, this 
is mainly due to the downfolding of the on-site e g character, as is explained in appendix 
A. The names and signs of the orbitals in subcells 2, 3, and 4 are thus the ones shown 
in the middle row of the figure. This is the natural choice when the structure is nearly 
cubic. 

3.1. Cubic t 2g bands 

In the cubic structure, the xy, yz, and xz Wannier orbitals are equivalent. Moreover, 
they are nearly independent because each orbital is even around its own plane and odd 
around the two other cubic planes. For that reason, the on-site and the lst-nearest- 
neighbour couplings between different orbitals vanish. As regards the coupling to 2nd- 
nearest neighbours, the yz orbital, for instance, does couple to xy orbitals at the four 
101 sites and to xz orbitals at the 110 sites, but the strength is negligible (6meV). 
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Figure 8. B^g xy, yz, and xz NMTOs for LaTiC>3. See also caption to figure 5. In 
the top row, the three different orbital shapes are in subcell f . In the bottom row, the 
mirror and the glide- mirror have been used to place them also in subcells 2, 3, and 4 
(see caption to figure 2). In the middle row, the orbitals have been grouped together 
in columns, and the signs have been chosen so as if the orbitals in subcells 2, 3, and 4 
were mere translations of the ones in subcell 1. We shall use the middle-row convention 
for naming the xy, yz, and xz orbitals in subcells 2, 3, and 4. This convention is the 
natural one when the distortion from the cubic structure is small. If the structure were 
cubic, the three orbital shapes would be identical, and the sum over the 4 orbitals in 
each of the middle-row pictures would give a Bloch wave with k=0. 

Before using the Wannier functions in figure 8, let us start simply by explaining 
the red t<i g band structure, e XiXj (k) , of cubic SrV03 in figure 3 in same way -but now 
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for each Bloch vector- as we discussed the gross features of the densities of states in 
figure 1; that is, in terms of a large basis set, the one giving rise to the black bands. 
In terms of Bloch waves, hybridization with oxygen states pushes t 2g band states up in 
energy and hybridization with A=Sr states pushes them down. In the cubic structure 
there is only one ABO3 unit per primitive cell and the Brillouin zone (BZ) is the blue 
one shown in figure 9. We may use the middle row in figure 8 to illustrate the three 
degenerate Bloch waves: with R being the positions of the B atoms (cells), each cell 
should be decorated with a phase factor expik • R; i.e. expik y in cell 2, expik z in cell 
3, expi (k y + k z ) in cell 4, etc. 

At T (k=0), there is no coupling to oxygen p because each t 2g Bloch wave is odd 
around the OA planes perpendicular to the plane of the wave and the O pn character 
is even. The direct ddn coupling is seen to be antibonding and it therefore tends to 
increase the energy at T, but since T marks the bottom of the band, the direct coupling 
cannot be the dominant source of the band dispersion. Finally, we realize that at T 
all A characters vanish, except the same Ad(t 2g ) character as that of the Bloch wave. 
Hence, e x . x . (0) is pushed down by interaction with Ad xx .. 

Going now from T to X 007T in figure 3, expik • R changes sign on sites 3 and 4, 
whereby the yz and xz waves become even around the horizontal AO mirror plane and 
can couple pdn antibonding along the z direction. This raises their band energy by 
^tpdn/ ( £ — e p) ~ 2 eV, minus the direct ddn contribution, and plus the lost contribution 
from the bonding interaction with respectively Ad yz and Ad xz , which is now forbidden 
and is not substituted by any other coupling to As or d. The xy wave becomes odd 
around the horizontal AO mirror plane, whereby its energy increases by merely a few 
tenths of an eV, due to direct dd5 hopping counteracted by the increased bonding 
interaction arising from the Ad xy character being substituted by Ad xz and Ad yz . 

Proceeding then from X 007T to M 7r07r in figure 3, the sign on sites 3 and 4 remains 
opposite to that on sites 1 and 2, but the sign alternates on the sheets which are 
perpendicular to x and not shown in figure 8. The energy of the yz band decreases by 
some tenths of an eV because the bonding interaction with Ad xy is now allowed, but 
this is counteracted by dd5. The energy of the xz band increases by less than 0.5 eV as 
the result of a further increase by ^ti dn / (e — e p ) , and of reductions due to oxygen pp 
hopping, ddir hopping, and coupling to A s. The energy of the xy band finally goes up 
by essentially 4^^/ (e — e p ) and becomes degenerate with the yz band. 

At R -Kirn the sign on all 6 nearest B neighbours is reversed from that of the middle 
row in figure 8, and the energy of the three degenerate bands is slightly higher than 
that of the xz band at M nOn. This is because of direct ddS hopping and because the 
coupling to A s is lost. Going finally from R back to T, the three bands stay degenerate. 

If we neglect the very weak inter-orbital coupling (6 meV in SrVOs), the dispersion 
of a cubic t 2g band can be written as a Fourier series: 

e m (k) = (k) = CI + 2 C, z m cos (xk x + yk y + zk z ) , (1) 

xyzj^xyz 

where x, y, and z run over all integers and where we have used the inversion symmetry of 
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Figure 9. Primitive cells of the monoclinic, nearly cubic and the quadrupled 
orthorhombic Pbnm structures (left), as well as their respective (blue and red) Brillouin 
zones (right). See also figure 2. Taking the pseudo-cubic lattice constant, |x| , as 
unity, the primitive reciprocal-lattice translations in the monoclinic structure are 
given by: 2ttG x = {[1 + (a + f3) /2] x + [{a - p) /2] y} / (1 + a + p + a(3) , 2nG y = 
{[(a - p) /2] x + [1 + (a + p) /2] y} / (1 + a + p + a/3) , and 2ttG z = z/ (1 + 7) . With 
k = k x G x + k y Gy + k z G z , the k x k y k z coordinates of the points marked in the blue, 
primitive monoclinic BZ are: G(r) 000, X 7r00, M 7T7tO, R 7T7T7r, and the equivalent ones. 
The high-symmetry points of the red, folded-in orthorhombic BZ are: G (r o ) 000, Z D 
OOf , X 0, Y |f 0, S ttOO, R 7r0§ , U | , T Q f § f , and the equivalent 
ones. The points explicitly listed are those used in figure 17. 

the orbitals and the lattice to combine the complex exponentials into cosines. In table 1 
we give the energy of the three degenerate t 2g Wannier functions, = ty° z ° yz = t°® z ° xz = 
e t2g —Spi and their transfer- or hopping integrals, such as the most important lst-nearest- 
neighbour integrals, tj^ xy = t^ z yz = t°^ xz = t w and t x Q y Xy = t y ° z ° yz = t° xz xz = t$, and 
the most important 2nd- nearest-neighbour integral, t xy ° xy = t yz ^ yz = t l xz xz = t' a . The 
description in terms of Wannier functions will be considered in detail in the following 
subsection, where we shall also explain why all the above-mentioned, most important 
hopping integrals are negative. Neglecting all other hopping integrals, the Fourier series 
(1) simplifies to: 

£ x . x . (k) = e t2g + 2t n (q + Cj) + 2t s c k + 4t' a CiCj, (2) 

where Xi = x, y, or z, and i, j, and k are all different. Moreover, q = cosk Xi . In this 
approximation the T— R bandwidth is —8t n — Ats, which for SrVC>3 is 2.4 eV. Hops to 
farther neighbours give another 0.5 eV. 

Had it not been for the direct and A-mediated effective ddS couplings between 
t 2g Wannier functions expressed by t$, the dispersion in figure 3 would have been two 
dimensional, with the bottom at T, top at M, and saddle-points at X. The latter would 
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to', to 
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give rise to a logarithmic van Hove singularity in the density of states. The strong 
reminiscence of this 2D singularity for SrVC>3 is clearly seen at the bottom left-hand 
side of figure 1. However, as we go along the series -but stay in the cubic structure- 
this sharp peak is gradually washed out due to increase of the A-mediated effective dd5 
coupling and, in the titanates, to A d bands overlapping the t 2g bands. 

3.2. The (xy,yz,xz) set of Wannier orbitals in the presence of GdFeO^-type distortion 

Distorting now the structures to the real ones, the density-of-states profiles are seen to 
sharpen in figure 1. This is because the increased O-A covalency pushes the Ad bands 
up above the t 2g band. But new structures arise, and these are caused by lifting the 
degeneracy between the xy, yz, and xz Bloch waves, by coupling between the waves, 
and by quadrupling the primitive cell. Let us now explain this -without getting into 
the details of the 4x3 orthorhombic bands- in terms of the (xy,yz,xz) set of Wannier 
functions for LaTiC>3 shown (before orthogonalization) in figure 8. 

In the top row of this figure we see that the t 2g orbitals have tc antibonding oxygen 
p character. This is the partner to the it bonding B t 2g character of the oxygen p orbitals 
in the cube faces seen in figures 5 and 6. Since those are the oxygen orbitals with which 
the Ad orbitals interact strongly (a interaction) and thereby cause the GdFeC>3-type 
distortion, they mix less with the B t 2g orbitals and the oxygen p character of the B t 2g 
orbitals in figure 8 is correspondingly weak. The oxygen character of the t 2g orbitals 
therefore decreases with increasing O-A covalency and the concomitant decrease of the 
shortest O-A distances (figure 7). Adding to this trend comes the reduction of oxygen 
character, cx ti d ^/ (e^ — e p ) 2 , caused by the fact that the 3d level of the earlier transition- 
metal ion, Ti, is higher above the O 2p level than that of the later transition-metal ion, 
V. This increase of td — e p seen in figure 1 is partly compensated by an increase of t p dn 
caused by the expansion of the 3d orbital when going from V to Ti. 

In addition to the oxygen pdn antibonding character, the t 2g orbitals have some 
A character and, in the distorted titanates, some e g character on neighbouring B sites 
where it bonds to the appropriate p lobe of a displaced oxygen. Such characters are 
barely seen in figure 8, but they are seen in figure 15 where a lower contour was chosen 
- albeit in subcell 2 and for a particular linear combination of t 2g orbitals. 

The oxygen character is decisive for the width of the t 2g band, that is, for the overall 
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size of the matrix elements of the Hamiltonian between t 2g Wannier orbitals on different 
sites (hopping integrals). The A character is important for lifting the cubic symmetry of 
on- and off-site matrix elements, and the B e g character induces on-site coupling between 
the orbitals. In order to demonstrate the specifics of this, we must first explain how the 
shapes of the Wannier orbitals are modified by the A spd and B e g characters: 

For all three t 2g orbitals in the top row of figure 8 (subcell 1), the red lobes 
attain some bonding sp and d 3z a i = (d xy + d yz + d xz ) j\/?> character on the nearest 
A neighbours, the ones at approximately \\\- The red lobes thereby stretch along the 
short A-B-A diagonal (orange bar in figure 7). This is not merely an effect of A-B 
covalency, but the oxygens set the stage as follows: the antibonding oxygen characters 
of the t 2g orbitals are deformed by bonding interaction with their nearest A neighbour 
(blue arrows in figure 7), and those A neighbours are thereby prevented from bonding 
to t 2g . As a result, the blue oxygen p lobes push the red t 2g lobes towards those two A 
ions which are not nearest neighbour to any oxygen, i.e. those along the [111] diagonal. 
The latter A ions are the ones free to bond with the red t 2g lobes, and they do this by 
contributing A sp d 3z 2^_ 1 character. This stretching of the t 2g orbitals along the shortest 
A-B-A diagonal is less pronounced in YTi0 3 than in LaTi0 3 , because in YTi0 3 , 01 is 
unusually close (red arrow in figure 7) to the A ion at \\\ and thereby partly blocks 
the A-B covalency. We shall later see how this difference between the two titanates 
influences their physical properties. 

The yz orbital is special in that its blue lobes attain bonding sp and d xy character 
on the two 2nd-nearest A neighbours, the ones at \-^\ in subcell 1. This is because the 
yz orbital has no p character on those oxygens, ^00, which have the \-^\ ions as their 
closest A neighbours, and therefore cannot prevent the blue yz lobes from hybridizing 
with them. At the same time, a blue yz lobe is pushed towards the ion at \-^\ by the 
red 02 p z lobe tending to bond with its nearest A neighbour, the one at \ = £\- The 
hybridization is with Ad xy because the tilt has moved the blue yz lobe towards the c 
direction, but the red 01 p y lobe then repels the blue yz lobe, which finally ends up 
running parallel to the flat Ol-A face (mirror plane). The strong displacement of 01 
leaves space for the blue yz lobe to bend over and become d xy -\ike at A \-^\- This effect 
is stronger in YTi0 3 than in LaTi0 3 , where also the 2nd-shortest A-B distance is 5% 
shorter than the average, rather than merely 1%. As a consequence, there is more Ad xy 
character in YTi0 3 than in LaTi0 3 . This difference between the two titanates will also 
turn out to be important for their physical properties. 

The xy orbital is the one whose orientation is most influenced by the GdFe0 3 -type 
distortion, simply because the axes of tilt and rotation are perpendicular to its lobe axes. 
The tilt moves the blue lobes towards the 2nd-nearest A ions at where the xy 

orbital attains some bonding A sp d xy character, but this is blocked by the attraction of 
the red 02 p y lobe to its closest A neighbour, as was explained above. For small GdFe0 3 - 
type distortions, i.e. in CaV0 3 , this reshaping of the t 2g orbital by 0-A covalency can 
be neglected, and the pure tilt of the xy orbital and its hybridization with A sp d xy are 
the most important effects of the distortion. 
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Finally, we restate that an oxygen p lobe which points opposite to the direction of its 
displacement in the cube face, may attain bonding e g character on the B neighbour. This 
enhances the p lobe and bends it outwards. This effect increases with the displacement 
of oxygen, i. e. with the degree of GdFe0 3 distortion. 



3.3. Effective hopping integrals 

In table 1 above, and in tables 2-5 we give the matrix elements of the LDA Hamiltonian 
(relatively to £p), 

LDA — / ± \njLDA I ± \ _ + xyz 



TTl^JJA — I _L \n_lLDA _ _L \ _ j.xyz 

n Ora',Rm — \ Xot 2g m' \ rL t F\ AR( 29 m/ — L m',m ' 

in a t 2g Wannier representation (A. 3). Here, R = x~R x + y~R y + zH z is a B site (see 
caption to figure 2) and x, y, and z are integers. The matrix element between orbitals 
m! and to, both on site 000, is t^,° m , and the hopping integral from orbital m! on site 
000 to orbital to on site xyz is t^, z m . Having the symmetry of d functions, these orbitals 
are even with respect to inversion in their own B site and, as a consequence, 

.xyz _ .xyz /n\ 
L m',m ~ L m',m- \° ) 

In this set of tables, the cubic basis (to = yz, xz, xy) is used, with orbitals named as in 
the middle row of figure 8. The values of the matrix elements, 

ttLDA / I \njLDA _ II \ _ ,R',R 

fl RW,Rra — \AR> t 2g m> \ n E F \ Xnt 2g m/ — l m',m 

with R' ^ 0, are given by the following rules, obtained by use of the mirror (z <-> —z) 
and the glide-mirror (x «-> y) : 

, (0,0,1), (x,y,z+l) _+xyz .(0,l,l),(x,»+l,z+l)_ ,yxz 

m' ,m ml ' (z^—z),m(z—+—z) m',m m'(x^y,z^—z),m(x^y, z-^—z) 

(4) 

,(0,0,0), (x,y,z) _ ,xyz .(0,l,0),(x,j/+l,z) _ ,yxz 

m',m m',m m' ,m m'(x^y), m(x+-*y) 

Before discussing these results, we stress that our downfolding to t<i a Wannier 
functions and evaluation of the Hamiltonian matrix is numerically exact. That is, 
our Wannier functions span the Kohn-Sham t 2g eigenfunctions exactly. Specifically, 
all partial waves other than B t 2g are downfolded, and their dependence on energy over 
the range of the t 2g band is properly represented. This is explained in appendix A and 
is illustrated in figures 3 and 10. The downfolding procedure used in model calculations 
such as [25, 29] is different, and not exact in the above-mentioned sense. In model 
calculations, only orbitals such as O 2p, and recently A 5d, expected to be relevant are 
downfolded, and their hopping integrals are assumed to follow the Slater-Koster rules. 
The latter may be a serious approximation, because, e.g., the O p orbitals have A and B 
s and p tails unless those orbitals are also being downfolded explicitly. Moreover, model 
calculations use Lowdin downfolding with e = e t2g (see equation (A. 5)). In the NMTO 
scheme, this corresponds to taking N — 0, that is, to using merely a single energy point. 
Once the Coulomb correlations are taken into account, not all the hopping integrals we 
tabulate may be relevant, of course, but the dominant ones are, and their values matter. 
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We now discuss our results and start by considering various contributions to a 
dominant t 2g hopping integral, namely the effective, lst-nearest-neighbour ddn hopping 
integral t w . A particularly important example will turn out to be t y )0 z 1 yz . This is the 
matrix element of the LDA Hamiltonian between the yz Wannier orbital in subcell 1 
and the similar orbital in subcell 3, defined as in the middle row of figure 8. Its value 
is seen to be negative (bonding), and to take the values — 281 me V in cubic SrV0 3 , 
-240 meV in CaV0 3 , -193 meV in LaTi0 3 , and a mere -65meV in YTi0 3 . This 
hopping integral between Wannier orbitals may be interpreted in terms of the energies 
of, and the hopping integrals between the simpler orbitals of a large basis -like the one 
giving rise to the black bands in figure 3 and like the ones discussed in for instance 
references [25, 29, 45]- containing Op, the Bd and As and Ad orbitals, all assumed 
to be orthonormal and to have interactions between lst-nearest neighbours only. If 
we perform Lowdin downfolding (A. 5) of the Hamiltonian in the large basis and use 
the middle row of figure 8 to keep track of the geometry and the signs, we obtain the 
following kind of result: 

+ 2 A+ 2 A+ 2 A+ 2 A+ 2 

, ,001 . °pdn ^ L Byz,As ^ Byz,Ayz ^ L Byz,Axz ^ Byz^Axy 

1 tt ■ l V z,yz ~ l ddiT 1 1 

y ' y e-t v e-e As e - e Ad e - e Ad e- e Ad 

i 2 At 2 At 2 At 2 At 2 

_ , pdit ^ L Byz,As ^ L Byz,Ayz ^ L Byz,Axz ^ L Byz,Axy /_>. 
— tddw 1 1 • ) 

£ — £ P €As — £ £Ad — £ tAd — £ £Ad ~ £ 
In the second line all signs are explicit because the direct hopping integral t dd7T is positive 
(anti-bonding), and e p < e < e As ~ e Ad when e is in the t 2 gband. The dominant term 
is the second one, the pdir hopping via oxygen. This term is bonding (negative) and 
is weakened by the direct ddn term, as well as by hops via A orbitals which are even 
with respect to the mirror plane between the two yz orbitals. The oxygen-mediated pdir 
hopping is strengthened by hops via A orbitals which are odd with respect to the mirror 
plane. 

That the first two terms in equation (5) have opposite signs, and that the second 
term dominates, should warn against interpreting hopping integrals between Wannier 
orbitals as overlap integrals. Moreover, characters barely visible in figures like 8 and 1, 
can contribute significantly to hopping if their energy is high, because the contribution 
to a hopping integral by a character of magnitude t 2 BA / (e A — e) 2 is t 2 BA / (e A — e) , that 
is, e A — e times larger. 

Now, the contributions from the four nearest A neighbours to an A term in 
expression (5) are identical only if the structure is cubic. In this case, actually, the two 
last terms in equation (5) cancel. Moreover, only if the structure is cubic do hopping 
integrals like t ddn and t pdv have their full value. Otherwise, they are reduced due to 
misalignment. 

In SrV0 3 the contribution from hops via A ions is small, as evidenced in figure 1 
by the weak Srd character in the Vt 2g band. But in the titanates, the t 2g band has 
almost as much Ad as Op character. This is not so evident from the appearance of a 
single t 2g orbital in figure 8, because for a Bt 2g orbital, B:A=1:8 rather than 1:1, as in 
the formula unit, AB0 3 . Similarly, for a Bt 2g orbital, B:0=l:4 rather than 1:3. In the 
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Table 2. H^ DA in meV for CaV0 3 [43]. 
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Table 3. #£ DA in mcV for LaTi0 3 using the older data [44]. 
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heavily distorted YTi03, the hopping integral tsyz,Axy to the 2nd-nearest Y, the one at 
position \-^-\ in subcell 1, is large, and that is the main reason for the anomalously small 
value of ty° z ^ yz in YTi0 3 . The values of the five other t n integrals, t x ° z 1 xz , t y l° yz = t l xzxz , 
and t x ° y ° = t° xy °, are normal and approximately —170 meV in this material. 

When estimating the value of t yz ° yz it should be remembered that, with the notation 
of the middle row of figure 8, the yz orbital in subcell 2 has the shape of the xz orbital 
in subcell 1; this is the reason why t yz ° yz does not have an abnormal Ad xy contribution, 
and also the reason why t yz ° yz = t l xzxz . 
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Table 4. H^ 2 DA in meV for LaTiOg [12].* 
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* We used their room-temperature data. For calculation of exchange constants, we also 
used the 8K data for LaTi0 3 . The result is shown in table 7. 



Table 5. H^ 2 DA in meV for YTi0 3 [20]. 
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The average of the t n integrals in the tables decreases by a surprisingly constant 
ratio of 22% for every step we proceed along the series SrV0 3 -CaV0 3 -LaTi0 3 -YTi0 3 . 
This holds for the recent LaTi0 3 data [12], while the older data [44] yield a 16% drop 
from CaV0 3 and a 29% drop to YTi0 3 . 

When deriving expression (5) from equation (A. 5) in configuration-space repre- 
sentation we have used low-order perturbation theory and have not properly inverted 
(P\e — Ti\P) . This means that the reduction of oxygen character in the Bt 2ff Wannier 
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orbitals due to 0-A covalency, and also the direct oxygen pp hopping, are not included 
in expression (5). But effectively, they both reduce the magnitude of all terms, except 
the first. Moreover, higher-order hopping processes such as B-O-A-B and B-O-A-O-B 
are neglected. Finally, the enhancement and outward bending of a p lobe on a displaced 
oxygen due to acquisition of bonding e g character on the B neighbour is also neglected, 
but it effectively increases t pdw towards that B neighbour. 

The second largest intra-orbital hopping integral is the effective 2nd-nearest- 
neighbour older integral, t' a . Like t n in (5), it is bonding (negative), but smaller because 
the distance is longer and because there is no oxygen-mediated hopping, except in case 
of the heavily distorted YTiOs. The typical integrals t xy ° xy and t xy ° xy are both —96 meV 
in cubic SrVC>3, —98 and —85 meV in CaVC>3, —47 and —46 meV in LaTiOs, and —50 
and +8 meV in YTiC>3. From the middle left part of figure 8, we gather that the integral 
txy°xy (txy°xy) f rom subcell 1 to the same subcell translated by b (a) has contributions 
from direct dda like hopping -which is now bonding and therefore strengthening,- from 
A mediated hopping, and from dpit-pp-pdit hops. For the two former, we get: 

9^2 oj-2 oj-2 oj-2 

/ _ no / ,110 \ , , Z ' L Bxy,As _ z ' L Bxy,Ayz _ ^ L Bxy,Axz ^ Bxy,Axy 

l a ■ l xy,xy \ l xy,xy ) ~ l dda + ~T 

v / £ — tAs £ — ^Ad £ — ^Ad ^ — tAd 

9^2 9J.2 9J.2 9J.2 

i, I ^ L Bxy,As ^ L Bxy,Ayz . ^ L Bxy,Axz z ' L Bxy,Axy /~\ 

— — \tdda\ 1 1 • W 

CAs — £ tAd — ^ tAd — ^ £Ad ~ ^ 

For t xy ° xy , the two important A ions are \\^- when seen from subcell 1, and ^^4^- 
when seen from b. For t xy ° xy , the two important A-ions are |-^-^r and ~y\^- as seen 
from respectively the origin and a. Keeping in mind the A characters of the xy orbital, 
we realize that although the intermediate A ions are different in cases b and a, strong 
A d xy or d yx / xz character does not exist on the same A ion; hence, A d mediated coupling 
is weak in all cases. The two first terms in expression (6) seem to dominate, except in 
YTi0 3 where the rotation of the xy orbital and the deformation of its d and p lobes 
are so strong that there is hopping, t pd , from the red p y lobe on 02 to the blue xy lobe 
translated by a. This gives a positive term, t p d n t p d/ (s — e P ) , which diminishes t xy ° xy . 

"xy,xy 



Due to the different deformation of the xy lobe, no such term exists for ^ 10 



The effective lst-nearest-neighbour <i<i<5-hopping, t$, e.g. t y ® yz , is an order of 



magnitude smaller than t w . It consists of the direct term — \tdd.s \ , plus A mediated 
contributions. It is negative for SrVC>3, changes sign near LaTiOs, and is positive for 
YTi0 3 . 

If the semi- quantitative expressions (5) and (6) for the effective hopping integrals 
are inserted in equation (1) for the dispersion of a cubic ti g band, the results of the 
discussion at the beginning of this section will of course be reproduced. 

We may define the root-mean-square value, t rms , of the hopping integrals by 
^t 2 = 12 (t rms ) 2 , where the sum runs over all neighbours and all three orbitals, and 
where 12 is the number of orbitals times the number of strong (t n ) hopping integrals 
in the cubic structure. These rms values are listed in table 6 and are seen to decrease 
with 19% from SrVC>3 to CaV03, with 15% from CaVC>3 to LaTiC>3 (using the recent 
data [12]), and by 11% from LaTiC>3 to YTiC>3. Even though the hopping integrals of 
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type ts and t' a , which decrease rapidly through the series, are included in the t 2 sum, the 
19-15-11% decrease of t rms is significantly smaller than the 22-22-22% decrease of the 
average t n . The reason is that the hopping between orbitals increases; without including 
inter-orbital hopping in the t 2 -sum, t rms would be 11% smaller in LaTi0 3 ([12]) and 
19% smaller in YTi0 3 . 

Table 6. Rms values of the hopping integrals in meV 

SrV0 3 [42] CaV0 3 [43] LaTi0 3 [44] LaTi0 3 [12] YTiO 3 [20] 
t rms 298 250 228 217 196 



3.4- Effective on-site matrix elements 

The Hamiltonian matrix elements most influenced by the lowering of the symmetry in 
the distorted perovskites are not the hopping integrals, but the on-site elements. The 
change of a hopping integral, however, usually perturbs the LDA bandstructure more 
than the same change of an on-site element, because the perturbation via the former 
is multiplied by the number of neighbours reached by that kind of hop. On the other 
hand, Coulomb correlations will effectively enhance the on-site elements and reduce the 
hopping integrals. 

As seen in table 2 for CaV0 3 , the energy of the xy orbital is 75 meV below that 
of the nearly degenerate yx and xz orbitals. One reason is that, for weak GdFe0 3 - 
distortion, the xy orbital is the one whose orientation is most optimized with respect to 
the A anions, specifically the two along [111]. The electrostatic field from these positive 
ions lowers the energy of orbitals pointing towards them. And so does the ligand field 
caused by hybridization with A orbitals, because for a diagonal matrix element, all 
interactions with characters of higher energy are bonding and push the energy of the 
Wannier orbital down. Conversely, all interactions with characters of lower energy are 
antibonding and push the energy up. As mentioned in connection with equation (5), 
tA > £ > e p . An additional reason for the energy of the xy orbital being the lowest is, 
therefore, that it has the least antibonding oxygen character. This comes about because 
the xy orbital interacts with four 02 ions, while each of the yz and xz orbitals interact 
with two 01 ions and two 02 ions. In addition, since Ca bonds more to 02 than to 
01 (the distance is 2% shorter), less 02 character is left for the V xy orbital. This is a 
pure ligand-field effect, since from oxygen there is no electrostatic effect due to the lack 
of JT distortion in CaV0 3 . 

In the titanates, the energy of the yz orbital at sites 000 and 001 is ~40 meV below 
that of the xy orbital. This is so because the orientation of the yz orbital towards 
the 1st- and 2nd-nearest A ions, those along the [111] and [111] diagonals, respectively, 
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Table 7. H^ DA in mcV for LaTiOg. 
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XV z 


000 


001 


010 


100 


000 


001 


010 


100 


000 


001 


010 


100 


m! , m 


























yz,yz 


465 


-195 


-185 


3 


420 


-68 


-168 


21 


471 


-41 


-211 


41 


Ob Z - Ob Z 


540 


-208 


3 


-185 


540 


-148 


21 


-168 


638 


-184 


41 


-211 


xy,xy 


500 


-21 


-185 


-185 


457 


-7 


-207 


-207 


566 


-10 


-167 


-167 


yz, xz 


-52 


-45 


76 


76 


36 


-89 


53 


53 


94 


-79 


73 


73 


xz, yz 


-52 


-45 


-40 


-40 


36 


-89 


5 


5 


94 


-79 


-35 


-35 


yz,xy 


-93 


-61 


-36 


-32 


-96 


-66 


-40 


-57 


-129 


-82 


-61 


-13 


xy,yz 


-93 


61 


-56 


73 


-96 


66 


-46 


95 


-129 


82 


-75 


58 


xz, xy 


-51 


50 


73 


-56 


-50 


82 


95 


-46 


-23 


91 


58 


-75 


xy, xz 


-51 


-50 


-32 


-36 


-50 


-82 


-57 


-40 


-23 


-91 


-13 


-61 



exploits both the electrostatic and the ligand fields. An additional ligand-field effect is 
that, in the titanates the bonding of the A ion to 01 is stronger than to 02. Therefore, 
less antibonding 01 character is available for the Ti yz and xz Wannier functions. This 
can actually be seen in the top row of figure 8. Finally, the B neighbour e g character 
which binds to the red back lobe of the 01 p y orbital -the one which is displaced most 
along its own direction- also contributes to lower the energy of the yz Wannier orbital. 
The orbital with the highest energy at sites 000 and 001 is xz. In LaTi0 3 , its energy is 
85 meV above that of the yz orbital, and in YTi0 3 , it is as much as 230 meV above. 
The main reason is that the xz orbital is the one of the three which is least favourably 
oriented with respect to the A ions, particularly in YTi0 3 . In addition, the xz orbital 
has little B e g character because its oxygen p lobes are not directed towards the face 
centres. A final reason why, particularly in YTi0 3 , the energy of the xz orbital is very 
far above that of the yz orbital is the JT distortion: the distance to Ti 000 of the 
oxygens along z and x is 3% shorter than the distance along y, and the pain character 
of the t 2 g functions is antibonding at the same time as the electrostatic field from the 
oxygen ions repel the near t 2g lobes. 

There is a strong on-site coupling between the two lowest orbitals, xy and yz, in 
the titanates. It is —89 meV in LaTi0 3 and —103 meV in YTi0 3 . This is due to the 
deformations of the xy and yz orbitals towards the 1st and 2nd-nearest A ions, those 
along [111] and [111]. With the signs chosen for the two orbitals, the lobes have the 
same sign when they point towards a near A ion, and the opposite sign when they point 
towards a far A ion. For this reason, both electrostatics and hybridization with all A 
characters enforce each other to make this on-site inter-orbital interaction strong and 
negative. 
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In all cases considered, the electrostatic and ligand-field effects work in the same 
direction. In the following, we shall therefore refer to anisotropy of the on-site 
Hamiltonian matrix for the ti g Wannier functions as a crystal-field effect, regardless 
of whether its origin is more electrostatic or covalent. 

3.5. Influence of the Jahn-Teller distortion in YTi0 3 

Although it is now generally recognized that the GdFeC>3-type distortion sets up a 
surprisingly large crystal field in the titanates, it is not agreed upon to what extent the 
JT distortion is important for the physical properties of YTiOs [21, 24, 29]. In table 7, 
we have included H^ 2 DA calculated for a hypothetical structure of YTi0 3 in which the 
octahedra are perfect, i.e., they lack the 3% JT-elongation of the Ti-02 bonds in the 
y(x) direction in subcells 1 and 3 (2 and 4), and the orthorhombic lattice constants are 
as in the real structure [20]. As expected, the energy of the xz orbital at site 000 is 
lowered, but remains 120 meV above that of the yz orbital at the same site. But the 
on-site couplings have not decreased, and the coupling between the xz and xy orbitals 
has even increased. Undoing the JT-distortion also influences the hopping integrals, but 
that does not make them more cubic: t n decreases between xz orbitals and increases 
between xy orbitals. With other plausible hypothetical structures without JT distortion, 
we obtained similar results. 

Recent high-pressure measurements [9] have shown that the 3% JT distortion 
essentially disappears in the pressure region between 9 and 14 GPa. At the same 
time, the Y ions move even further away from their cubic positions. At 16 GPa, the 
Ti-02 distances in the x and y directions have become nearly equal, with the Ti-Ol 
distance being just 1% smaller. Moreover the GdFe03-type tilt has increased from ~20 
to 21°, while the rotation is unchanged. The on-site Hamiltonian and nearest-neighbour 
hopping integrals calculated for this structure are included in table 7. We see that the 
ordering of the orbital energies is the same as for the JT distorted and undistorted 
structures at normal pressure, and that the energy of the xz orbital is "merely" 167 
meV above that of the yz orbital, while that of the xy orbital is 95 meV above. These 
results are similar to the ones obtained for the hypothetical GPa structure without JT 
distortion, and so is the result that the on-site coupling between the yz and xy orbitals 
increases upon removal of the JT distortion and readjustment of the Y positions. 

In conclusion, the JT distortion influences the crystal field in YTiOs, but is not its 
source. We shall return to this subject at the end of section 3.6, as well as in sections 
3.7, 5.4, and 6. 

3.6. Influence of GdFeO^-type and JT distortions on the t2 g bandstructures 

In the preceding subsections, we have clarified how the on-site and hopping matrix 
elements of the LDA Hamiltonian develop along the series. Now we need to understand 
how the GdFe0 3 -type distortion perturbs the spatial coherence of the simple, cubic 
Bloch waves. 
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Table 8. t 2g edge-to-edge (W t2g ) and rms (W) bandwidths in eV. 




W 



SrV0 3 [42] CaV0 3 [43] LaTi0 3 [44] LaTi0 3 [12] YTiO 3 [20] 
2.85 2.45 2.09 1.92 2.05 

2.85 2.39 2.18 2.08 1.87 



The LDA bandstructures in the region of the t2 g bands are shown in figure 10 for all 
four materials. Like in figure 3, the red bands, £j (k) , are obtained from the t2 9 Wannier 
functions, e.g. by forming Bloch sums (see equation (A. 6)) of the Wannier orbitals and 
diagonalizing the Bloch transformed Hamiltonian, 



We thus obtain the eigenvalues, E{ (k) , and eigenvectors, Uj (k) . Since the orthorhombic 
structure has four (equivalent) subcells per translational cell, there are 4x3=12 
orthorhombic t 2g bands. In equation (7), T are the orthorhombic translations, and 
R and R' run over the four B sites, 000, 010, 001, and 011. The structure is specified 
in the caption to figure 2. The Hamiltonian is a periodic function of k in the lattice 
reciprocal to the orthorhombic lattice. The corresponding BZ is shown in red in figure 
9, whose caption specifies the reciprocal-space structure. The 12 bands are doubly 
degenerate on the faces of the orthorhombic BZ due to the presence of the mirror and 
the glide mirror in the space group. 

As we proceed along the series, the t 2g edge-to-edge bandwidth, W t2g , is clearly seen 
to decrease. It is tabulated in the top row of table 8. In the second row, we give an rms 
bandwidth, W, which we define to be proportional to t rms (table 6) with the prefactor 
chosen such that W = W t2g in case of cubic SrV0 3 , i.e.: W = (W t2g /t rms ) SrV03 x t rms . 
We see that the decrease of the edge-to-edge bandwidth does not follow the trend of 
t rms . So, clearly, band shapes change along the series. 

In view of the extreme simplicity of the cubic bandstructure of SrV0 3 (subsection 
3.1), the orthorhombic ones are bewilderingly complicated. They do not immediately 
tell us what the progressing GdFe0 3 -type distortion actually does. Simple theories have 
assumed that all hopping integrals scale uniformly with the bandwidth, i.e., that the 
bandstructures of all four materials look the same; clearly, they do not. 

The next simplest case is that of a pure, cooperative JT distortion, e.g. the one in 
which the energy of the yz orbital on site 000 is A lower than that of the xz orbital on 
the same site. The glide mirror, which exchanges yz and xz on the nearest neighbours 
in the x and y directions, will then make a yz/xz orbital feel a spatially alternating 
potential, =f|A exp in (x + y) = =f|A exp[i (ir, ir, 0) • R]. This makes the cubic e yz (k) 
and E yz (k— 7T7rO) bands split by A at the surface in k space where they would have 
crossed had there been no JT distortion. The same holds for the xz band. In this 




(7) 



T 
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Figure 10. Low-energy LDA bandstructures of SrV0 3 , CaV0 3 , LaTi0 3 , and YTi0 3 . 
The cubic bands for SrV0 3 (figure 3) have been folded into the orthorhombic BZ 
(figure 9). The black bands were obtained with a large NMTO basis, while the 
coloured ones were obtained with the truly-minimal B d fog) basis. The two are 
indistinguishable, except for LaTi0 3 where the bottom of the La 5d-band overlaps the 
top of the Ti ^g-band because the older structural data [44] was used here. With the 
recent, slightly more distorted structure [12], O-A covalency removes this overlap as 
discussed in section 2 and figure 1, and as shown in figures 4, 12, and 17. In order to 
emphasize the development of the bandstructure along the series, the gap between the 
lower 1/3 and the upper 2/3 of the bands has been shaded. 



case, the xy band is unaffected, and there is no need for folding the bands along the z 
direction, i.e. by the reciprocal lattice vectors 007r and tttttt. 

For the t\ g materials, however, the dominant crystal-field splitting, A = ty° z ° yz — 
t'xzxzi * s & t m ost a tenth of the bandwidth, and yet, individual bands are seen to shift by 
up to half the bandwidth when we proceed along the series in figure 10. The following 
shifts (splittings) are large and important, because they move bands which are near the 
Fermi level: the 4th band at T moves up, and in YTiC>3 it almost empties. At the same 
time, the 5th band at T and the degenerate 3rd and 4th bands at S Q move down, so 
that in the Stoner model for ferromagnetism in YTiC>3 they become occupied (see also 
figure 18). This downward shift of the 3rd and 4th bands near S D is the cause for the 
gradual development of a pseudo-gap between the lower 1/3 and the upper 2/3 of the 
t2g band. In comparison, the A-gapping of the yz and xz bands along the lines Y T 
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and X U ^| k z is less spectacular. 

Coulomb correlations will of course change all of this, but in a fairly simple way, 
at least in the DMFT as we shall see in section 5.3: the one-electron excitations can be 
described as the LDA Hamiltonian plus a self-energy, which in the single-site DMFT is 
simply an energy- dependent, complex on-site matrix, S mm / (u) . It is therefore important 
to understand the LDA bandstructure and, even better, to have an analytical model of 
it, to which a given £ mm / (u) can be added or even determined, say, by a model DMFT 
calculation. Such an analytical model should for instance be useful for understanding 
why recent de Haas-van Alphen (dHvA) measurements for CaV0 3 gave a Fermi surface 
more cubic than the orthorhombic one obtained from the LDA [46]. 

An analytical bandstructure model may be obtained by transforming from the 
mixed (k, R) representation (7), where R are the four sites in the orthorhombic cell, to 
a pure wave-vector k — Q representation, where Q are the four smallest, inequivalent 
vectors, 

Q = 000, ttttO, OOtt, and tttttt, (8) 

of the orthorhombic reciprocal lattice (Q ■ T = 2n x integer) . In figure 9, these are 
the blue points T, M xy , X 2 , and R of high symmetry in the cubic (actually primitive 
monoclinic) BZ. The transformed Bloch-waves are: 

4= Yl i R > m > k > e " Q R = E ~r= E l R ' m > e*( k - Q H R+T > 



R R 

1 



^|R,m)e i(k - Q) - R = |m,k-Q), (9) 

v ^ c R 

where, in the first expressions, R runs over the 4 sites and, in the last expression, over 
all (L c — > oo) translations of the cubic lattice. T runs over all (L Q ) translations of the 
orthorhombic lattice (see figure 2). |R, m) denotes a yz, xz, or xy Wannier orbital 
with the convention of the middle row in figure 8, so that the energy and shape of 
the m orbital depends on which of the four types of sites it is centred on. In the Q 
representation, the expression for the Hamiltonian matrix becomes: 

<ro', k - Q' | U\ m, k - Q) = £ e^> R , (10) 



R 



with 



i 4 

,Q,R _ 4 ,R',R'+R iQ R' 

L m',m — a / j i m',m ' ' 



4 

R' 



and where for simplicity of notation we have dropped the superscript LDA on H. Note 
that the difference between any of the four Q vectors also belongs to set (8). If t R 



m ,m 



were independent of R', as would be the case if the energy and shape of each m orbital 
were independent of its position, we would have: t^f^ = 5q t R m , and the Hamiltonian 
would have cubic translational symmetry. But the GdFe03-type distortion of the on- 
site energies and shapes of the Wannier orbitals introduces coupling between the four 
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orthorhombic subdivisions of the cubic BZ. Each subdivision is denoted by a specific Q 
vector in the set (8), or one of their equivalents. 

We must now calculate tz) m using the unitary matrix, 



3 iQR 



41 



X 



Q 


R 000 


010 


001 


011 


000 


1 


1 


1 


1 


7T7T0 


1 


-1 


1 


-1 


OOtt 


1 


1 


-1 


-1 


1T1T1T 


1 


-1 


-1 


1 



and the space group (Pbnm) and orbital symmetries expressed by equations (4) and 
(3). In order to obtain simple, explicit expressions, let us limit the lattice sum in 
(10) to include merely the on-site, the six lst-nearest, and the six t^-type 2nd-nearest- 
neighbour hoppings. For the latter, we neglect the tiny differences between the three 
numbers, t yzyz = t xzxz , t xzxz = t yzyz and ^ (t x 
of the Hamiltonian can then be expressed as: 



n,..,;, + tl x 1 y°xy) ■ The Q'Q matrix element 



(k-Q'|W|k-Q) = t 



Q'-Q,000 



-2t 



Q'-Q,001 



-2t 



Q'-Q,010 



C'ii 



■2f 



Q'-Q,100 



C, 



+ 4f Q'^ 1 110^^ + 4i Q'^,101^^ + 4t Q'-Q 1 011 Cw ^ + 4i -Q'-<J ) 110 Sa!Si/> (U) 

where c x = cos (k x — Q x ) and s x = sin (k x — Q x ) , etc. In order to specify the t Q ' R 
matrices in orbital space, it will finally prove convenient to define the following unit, 
symmetric, and anti-symmetric matrices: 



E yz = 





yz 


xz 


xy 


yz 


1 








xz 











xy 












Cjyz,xz = 





yz 


xz 


xy 


yz 





1 





xz 


1 








xy 












j^yz,xz = J^xz,yz 



yz 


xz 


xy 


yz 


1 





xz —1 








xy 









(12) 



etc. 



Note that the A matrices are defined to have 1 or in their upper triangle, and —1 or 
in their lower triangle. 

The Q = Q' = element, (k \7i\ k) , is a "cubically" averaged 3x3 Hamiltonian. 
This is a generalization of the non-interacting cubic bands given by (2). Its on-site and 
001-hopping matrices are given by: 



Q=000,R=00z 



Clz) (E yz + & 



I i00z fpxy , J.UU2 qyz,xz 
' xy,xy ' yz,xz 



(13) 



if z takes respectively the value and 1. Here, the coefficients are linear combinations 
of the basic Hamiltonian matrix elements given in tables 1-5. For instance is 

/M0z \ = I (fOOz , ,00z \ 
\ L LZ,Lz/ — 9 \ L yz,yz ' i xz,xz) 
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Figure 11. Cubically averaged (Q = 0) ti g LDA bandstructures in eV, in the primitive 
mono-clinic BZ (left), and folded into the orthorhombic BZ (right). The letters above 
the k x k y k z coordinates denote the orthorhombic high-symmetry points (red in figure 
9), k + Q, into which the corresponding cubic point is folded on the right-hand side. 
Bands with predominant xy, yz, and zx character are respectively red, green, and blue. 
Black bands are of strongly mixed yz and zx character. The definitions of characters 
are those of the middle row in figure 8. The bands were obtained with the parameters 
for Q=000 in table 9. In order to obtain the bandstructures in figure 10, albeit in the 
2nd-nearest neighbour approximation, the folded-in bands should be coupled using the 
matrix elements (14)-(16). 
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Table 9. Coefficients of the on-site energies and hopping integrals in meV (see 
equations (11)— (16)) 

R 000 001 010/100 

Sr Ca La Y Sr Ca La Y Sr Ca La Y 

Q 

000 (t L2 , L2 ) 625 616 488 490 -281-245-201-122 t UZ;llz -281-223-185-184 

t xy ,xy 625 542 486 417 -33 -23 -22 -8 t xy , xy -281-231-183-162 

tJ tLje -33 -17 1 28 

t yZ:XZ 4 -42 43 -21 -42 -63 (t nz ^) 6 16 15 

vrvrO t yZ: y Z -t xz , xz 8 -85-230 9 15 113 t yz , xz -t xz , yz 78 118 111 

00vr t yz , xy +t XZ:Xy 0-21-131-120 -11 -7 19 t nz , xy +t XV:llz -51 -88-118 

t^z,xy~^~txy,^z 8 40 43 

vrvr^ t yZ:Xy -t XZ:Xy 0-13 -47 -86 -61-111-147 t y l° xy -t xy ° yz 17 20 10 

tfz° xy -tfy% 58 106 87 

R 110/011/101 

Sr Ca La Y 

Q 

000 t' c -96-91 -47 -27 

s' 7 1 29 



the on-site energy or 001 t n hopping integral, averaged over the yz and xz orbitals. The 
matrices for hopping parallel to the mirror plane (R = 010 or 100) have the very similar 
form: 



IZ.L2 



y.Q=000,R=010/100 _ j. puiz _i_ j. pLZ , / rpxy 1 U \ C 

1 t iiz,iiz- C/ ~r i-l2,lz- C/ ~r l xy,xy I - J t \^hz,lz/ 

where 11 and l, respectively, denote the directions parallel and perpendicular to the 
direction of the hop, i.e., for R = 010, 11= y and i_= x, whereas for R = 100, 11= x and 
i_= y. So the first three terms are the t n or t$ hoppings in the plane, and the last term 
is the average of the hopping between the yz and xz orbitals along 010, or 100: 

It ) = -(t +t ) = - (t mo + t 010 ) = - (t 100 + t 100 ) 

\^i\z,^z/ — 2 \ v uz,uz ~ >>uz,iizj 2 \ yz,xz ' b xz,yzJ \ xz,yz ' L yz,xz) ■ 

Remember that = 7^ t°xz,yz = ^x°z°yz- The hopping to the 2nd- nearest 

neighbours will be diagonal in both orbital and Q spaces: 

t Q ' 110 = 5 Qi0 t' a E xy , t^ wl = 5 Qi0 t' a E™, t Q ' 011 = 5 Qi0 t' a W, t^ 110 = 5 Qi0 s' a E*», 
as a consequence of the approximation made above. Here, 



,110 ,,110 +o(t 011 + / 101 1 and <t' = - ( t 110 - t nu 

xy,xy ~ L xy,xy ~ \ L yz,yz ~ L xz,xzJ 5 a a r> \ xy,xy 



2 



"xy,xy 
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Note that in these expressions, and in (14)— (16) below, the t QR , E, S, and A matrices 
in orbital space are the only such matrices; the prefactors, such as t° xyxz , are parameters 
and their subscripts do not label elements of 3 x 3 matrices. 

The values of these parameters are displayed in table 9, whose top four rows, 
denoted Q=000, determine the cubically averaged bands. We see that in SrV0 3 and 
LaTi0 3 the average (on-site) energy of the xy band is the same as the common, average 
energy of the yz and xz bands, whereas in CaV0 3 and YTi0 3 , it is ~75 meV lower. The 
values of the hopping integrals confirm that the cubically averaged bands narrow along 
the series (due to misalignment and increased theft of oxygen character), and that the 
hopping in the z direction is anomalously small in YTi0 3 , which therefore has the most 
anisotropic t 2g band. The reason for that is the smallness of t y ° zyz due to hybridization 
with Y 4d xy , as was discussed in connection with equation (5). 

The cubically averaged bandstructures are shown in figure 11. We clearly see how 
these bands develop from the cubic, non-interacting, nearly two-dimensional xy, yz, and 
xz bands of SrV0 3 into three distorted, monoclinic bands in which nearly degenerate yz 
and xz levels split by ± [t y ° z ° xz + ^ yzxz cos k z + 2 (t„ z ^ z ) (cos k x + cos k y )] . The xy band 
(red) stays pure and nearly two-dimensional. For YTi0 3 , the weak k z dispersion of the 
yz band (green) is clearly seen in the (010) plane where there is no mixing with xz. 
The xz band (blue) will have the same behaviour in the (100) plane in this cubically 
averaged bandstructure. The letters r o , T G , S G , R Q , Z Q , and X Q on the abscissa refer 
to the orthorhombic high-symmetry points, k + Q (red in figure 9), into which the 
monoclinic k points (blue in figure 9) should be folded prior to Q-coupling. The folded- 
in bands are shown on the right-hand side, and comparison with the real ^ bands 
in figure 10 reveals that, in fact, the cubically averaged bandstructures go "halfway" 
towards the real bandstructures. 

The hybridization between two cubically averaged bandstructures displaced in k 
space by the orthorhombic reciprocal lattice vector Q=7r7r0 is given by the matrix (11) 
with: 

2t Q=^0,R=00, = ( f 00 Zyz _ t 00y {E y Z _ ^ {u) 

9-Q=7T7rO,R=010/100 _ (• _ • \ Ayz,xz 

^ \ L yz,xz L xz,yz) , 

and is seen to express the asymmetry between the yz and xz orbitals, as well as the 
coupling between them. Here, the top line provides the A splitting by (t y ° z ° yz — t x ° z ° xz ) + 
2 (ty°z]yz ~ tx°z]xz) c z of the s yz (k) and e yz (k— ititO) bands at their crossings. The same 
holds for the splitting of the e xz (k) and e xz (k— irirO) bands. Note that since t Q yzyz — t x Q zxz 
is a parameter, it should not be substituted by t x ° z z xz — t yz z yz when describing the xz 
band, because that change of sign is accounted for as a prefactor to E xz . Note also that 
y.010 _ +010 _ +100 _ +100 because t 010 = t 100 t 100 = t 010 

yz,xz xz,yz yz,xz xz,yz UCL-dUoe l y Z ^xz yz,xz I xz,yz xz,yz' 

Similarly, the hybridization between two cubically averaged bandstructures 
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displaced by Q=007r is specified by: 

rwQ=007r,R=000 _ / + 000 i ,000 \ / nyz,xy , nxz,xy\ 
Zl — \ l yz,xy "+" l xz,xy) l D + D I 

2f Q=00*,R=001 = (yOOl^ + yO<H J ( A yz,xy + ^sj,) (15) 

9-Q=00tt,R=010/100 _ ( f i y \ onz.xj/ , /y _|_ y \ cjLZ,xj/ 

^ v \ L \\z,xy ~ v xy,\\z) u ~ \ L i-Z,xy ~ L xy,i_z ) u i 

and is seen to express the average of the coupling from the xy orbital to the yz and 
xz orbitals. As above, n and l denote the direction parallel and perpendicular to the 
direction of the hop, i.e., for R = 010, n= y and i_= x, whereas for R = 100, n= x 
and L = y. Note that t^ xy = -f™,, that t™ xy = t™ xz ± t\% y = t™ yz) and that 

y010 _ ylOO / ylOO _ yOlO 

xz,xy ^xy,yz I ^yz,xy ^xy,xz' 

Finally, two cubic bands displaced by Q=7T7r7r hybridize by means of a matrix (11) 
specified by: 

2^.Q=7r7r7r,R=000 _ ^000 _ yflOO \ / gyz,xy _ gxz,xy\ 
\ yz,xy xz,xyj V / 

2 yQ™,R=001 = ( y001^ _ y001J {A yz,xy _ A xz,x V) (lg) 

9yQ=vr7r7r,R=010 _ /y010 _ yOlO \ Ayz,xy , /.010 _ yOlO \ Axz,xy 
™ V yz,xy ^xy,yz) ' \ L xz,xy L xy,xz) 

9yQ=7r7r7r,R=100 _ _ /.010 _ yOlO \ \yz,xy _ /y010 _ y010 \ Axz,xy 
" VxZjXj/ ^xy,xz) \ yz,xy ^xy,yz) 

(7 100 _ f 100 ^ 42/^2/ 1 C/100 _ ylOO \ Axz,xy 
Vyz,xy L xy,yz) ' \ xz,xy xy,xz) 

This matrix expresses the difference between the couplings from xy to yz and to xz. We 
emphasize that there is no Q coupling between xy bands. 

The values of the coefficients given in table 9 exhibit the general trend that the 
non-cubic perturbations increase along the series. We also see that details differ: which 
ones of the many non-cubic couplings dominate depends on the material. Moreover, the 
on-site, k-averaged splittings, t y Q zyz — t x ° zxz) and inter-orbital couplings, t y ® z ° xy ±t xz ° xy , are 
not much larger than their modulations given by the corresponding hopping integrals. 
This is very different from the situation in materials with strong JT distortions, where 
Q coupling due to on-site terms dominate. So according to the present conventional 
description of the ti g perovskites in terms of Q couplings, the distortion from the 
cubically averaged to the real bandstructure is the sum of several terms with varying 
signs. 

As a most important example, let us discuss the pseudo-yap. More specifically, let 
us estimate the levels at S c 7r00 with intermediate energy by using expressions (11)— (16), 
table 9, and figure 11: 

The relevant cubically averaged yz levels are those at k — Q = 7r00 — 007r = tiOti 
and 7r00 — tttiO = OnO. Their energies are: 

Syz (vrOvr) / e yz (OvrO) = <^° L2 > - < =p 2 - t nz , uz + t L2 , L2 ) 

= ( 1009 980 676 598 ) T ( -66 -58 -30 180 ) 

= ( 1075 1038 706 418 1 /f 943 922 646 778 1 meV 



for SrVC>3 through YTiC>3. Note that in YTiC>3, the e yz (7i07i) energy is anomalously low, 
as is also seen in figure 11. This is one reason for the pronounced pseudo-gap in that 
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material. The yz levels are separated by Q' — Q=7T7T7f and, according to (16), or table 
9, do therefore not couple. Degenerate with these yz levels are the cubically averaged 
xz levels, 

£ xz (Otttt) /e xz (ttOO) = e yz (ttOtt) /e yz (OttO) , 

at k — Q = 7r00 — tttttt = Otttt and 7r00 — 000 = 7r00. Also these xz levels are separated 
by tttttt, and do therefore not couple with each other. The yz level at ttOtt and the xz 
level at 07f7f are separated by ttitO and therefore couple via a matrix element obtained 
from equation (14), which is seen to vanish, however, because c x + c y = 0. The same 
holds for the yz level at 07f0 and the xz level at 7r00. Finally, the yz level at ttOtt and 
the xz level at 7r00 are separated by 00-7T so that, according to (15) or table 9, they do 
not couple. In conclusion, the four yz and xz states at S G with intermediate energy can 
only couple via xy states, which we consider next: 

When folded into S , all four cubically averaged xy bands (the red ones in figure 11) 
are nearly degenerate. One pair of degenerate levels come from k — Q' = 7r00 — 007r = 
ttOtt and 7r00 — tttttt = Otttt, and the other from k — Q' = 7r00 — 000 = 7r00 and 
7r00 — ttttO = OttO. As for the yz and xz levels, the two xy energies are those of nOn and 
7r00 points, specifically: 

[e xy (TT0TT)=e xy (Otttt)] / [e xy (ttOO) =e xy (OttO)] = t x ° y ° xy — 4t' a =F 2t 001 

-( 

= ( 1075 952 718 541 ) 943 860 630 509 ) 

Since there is no Q-coupling between xy bands, we now only have to couple xy to 

yz and xz. 

The yz level at ttOtt couples to the lowest, degenerate xy level, that is the level at 
7r00 and 0-7f0, because for those, Q' — Q = ttOtt — 7r00 = OOff and ttOtt — OttO = tttttt, 
respectively. It does, however, not couple to the xy level at ttOit, because for those states, 
Q' — Q = ttOtt — ttOtt = 000 and ttOtt — Otttt = ttttO. Similarly, the xz level at 07f7f, 
which is degenerate with the yz level at 7fOTt, also couples to the xy level at ttOO and 
0-7f0, because, here, Q' — Q = OtTtT — 7r00 = tttttt and OtTtt — 07f0 = 00tT. As for the yz 
level, the xz level does not couple to the xy level at Otttt. 

The yz level at OttO, which is lower than the one at ttOtt, except in YTiC>3, 
couples to the highest, degenerate xy level, the one at ttOtt and 07T7f, because for those, 
Q' — Q = OttO — ttOtt = tttttt and 07f0 — 07f7f = 007T. It does not couple to the xy level 
at 07f0. Similarly, the xz level at 7r00 couples to the xy level at 7rO-7f and 07f7f, because 
Q' — Q = 7r00 — 7r07f = 007r and 7r00 — 07f7f = tttttt. Again, there is no coupling between 
xz and xy levels at 7r00. 

All couplings relevant for the pseudo-gap at S G thus have Q' — Q = 007r and tttttt, 
and are therefore described by the five hopping parameters listed in the corresponding 
rows of table 9. We emphasize that the dominant crystal-field coupling, caused by 
tyTvz ~ t°xz°xz-> nas Q' — Q = ttttO, and is therefore not relevant for the pseudo-gap at S G . 



xy,xy i — xy,xy 

1009 906 67 J 525 ! i -66 -46 -44 -16 ) 

943 860 630 509 meV. 
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For the couplings via Q' — Q = 007r, we find from equations (11) and (15): 

(xy, nOO \H\ nOn, yz) = (xy, OnO \Ti\ Oim, xz) 

- I f+000 , t 000 \ , (,001 , .001 \ +(t +t \_ (t +t \ 

2 IVffl xz,xyj ~ Vyz,xy ~ h xz,xy) ~ \ b uz,xy ~ L xy,uz) \^LZ,xy ~ b xy,i_z) 

= i (-21-131-120) + (-11 -7 19) + (-51 -88 -118) - (8 40 43 ) 

= (-81 -201 -202) meV, 

for CaVC>3 through YTi03. These couplings are numerically large for the titanates 
because terms add up. The other 007r couplings are: 

(xy, Onn \TC\ 0n0, yz) = (xy, nOn \Ti\ nOO, xz) 

= I (t 000 + f 000 ) - (t 001 + t 001 )-(t +t ) + (t +t ) 

2 \ L yz,xy ' L xz,xy) \ L yz,xy ' L xz,xyJ \ l \\z,xy \ h xy,\\z) \ \ L LZ,xy i L xy,Lz) 

= i (-21 -131 -120) - (-11 -7 19) - (-51 -88 -118) + (8 40 43 ) 
= ( 60 70 82 ) meV, 

which are small for the titanates due to cancellation of terms. For the couplings via 
Q' — Q = nnn we find from expressions (16): 

(xy, 0n0 \Ti\ nOn, yz) = — (xy, nOO \7i\ Onii, xz) 

- i (t 000 — t 000 ) 4. (t 001 — t 001 ) — (t 010 — t 010 ) — (t 010 — / 01 ° ^ 

2 V yz,xy ^xz,xy) ' \^yz,xy ^xz,xy) \ yz,xy ^xy,yz) \^xz,xy ^xy,xz) 

= ^ (-13 -47 -86) + (-61 -111 -147) - (l7 20 lo) - (58 106 87) 

= (-143 -261 -287) meV, 

which are large for all three materials, and 

(xy, 7r07r \H\ 07r0, yz) = — (xy, Oirir \H\ n00, xz) 

_ 1 (.000 _ ,000 \ _ (.001 _ ,001 \ 1 1 .010 _ ,010 \ 1 / ,010 _ ,010 \ 

2 V y z i x V xz,xy) \ yz,xy xz,xy) ' \ yz,xy xy,yz) ' \ xz,xy ^xy,xz) 

= i (-13 -47 -86) - (-61 -111 -147) + (l7 20 lo) + (58 106 87) 

= (l30 214 20l) meV, 

which are large as well. 

Finally, we can form the matrix (7r00— Q' \H\ 7r00— Q) for the 8 states of 
intermediate energy at S G . This matrix blocks into two. The first block involves the 
states: \xy, 7r00— 000) , \xy, ttOO— ttttO) , \yz, 7r00— 007r) , and \xz, nOO— nnir) . Along its 
diagonal, it has the degenerate levels e xy (7r00) and e xy (07r0) , followed by the degenerate 
levels e yz (ir0TT) and e xz (0inr). The off-diagonal elements are: (xy, 7r00 \H\ 7r07r, yz) = 
(xy, 07r0 \H\ 07T7T, xz) and (xy, OnO \H\ irOn, yz) = — (xy, n00 \H\ Onii, xz) . The states 
contributing to the second block are: \zx, 7r00— 000) , \yz, nOO— irnO) , \xy, n00— 00n) , 
and \xy, 7r00— tttviv) . Along its diagonal, there are the degenerate levels e xz (ttOO) and 
e yz (07r0) , followed by the degenerate e xy (nOn) and e xy (Otttt) levels. Its off-diagonal 
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elements are: (xy, Otttt \H\ OttO, yz) = (xy, ttOtt \H\ 7r00, xz) and (xy, ttOtt \H\ OnO, yz) = 
— (xy, Otttt \7i\ ttOO, xz) . Each of these 4x4 blocks can now be transformed into two 
degenerate 2x2 blocks. For the first set of states, this 2x2 block is simply: 

e xy (ttOO) = e xy (OttO) \J (xy, ttOO \H\ 7r07r, yz) 2 + (xy, OttO \H\ ttOtt, yz) 2 \ 
he s yz (ttOtt) = e xz (Otttt) 



860 V81 2 + 143 2 \ / 630 V201 2 + 261 2 \ / 509 V202 2 + 287 2 
he 1038 ) ' I he 706 ) ' I he 418 



meV, 



where we have inserted the values for CaVC>3 through YTiOs. Note that the size of the 
hybridization is \/2 larger than a typical matrix element. Diagonalization yields the 
corresponding S levels: 

£l(So)= (me)' (1000)' (no') meV - 

The 2x2 matrix and eigenvalues for the second set of states are respectively: 



e xy (ttOtt) = e xy (Otttt) y (xy, Otttt \H\ OttO, yz) 2 + (xy, ttOtt \H\ OttO, yz)' 
he e yz (OttO) = e xz (tt00) 



'952 v / 60 2 + 130 2 \ / 718 v / 70 2 + 214 2 \ ( 541 V82 2 + 201 2 
he 922 J ' I he 646 J ' 1 he 778 



meV 



and 




We realize that these eigenvalues give the trends well, although accurate agreement with 
the bands in figure 10 requires inclusion also of the remaining 4 states at S G , as well 
as inclusion of longer-range hoppings. The decrease through the series of the lowest, 
degenerate eigenvalue, that is, of the 3rd and 4th bands at S D , from 943 meV in SrVC>3 
and 762 meV in CaVC>3 to 336 meV in LaTiC>3 and, finally, to 110 meV in YTiC>3 
is spectacular. Note also the inversion of the cubically averaged e yz (7r07r) =e xz (Otttt) 
and e xy (tt00) =e xy (OttO) levels in YTiC>3. In both titanates, the couplings between 
the xy and the yz and xz Bloch waves are strong, but the anomalously low-lying 
e yz (ttOtt) =e xz (tt00) in YTi0 3 makes the difference. 

The chemical reason for the development of the pseudo-gap is the increasing residual 
A character in the ti g band. This was mentioned before and is demonstrated for the 
titanates in figure 12, where we have projected the bands onto the Ad 3z 2^_ l and A 
d xy partial waves. The lowering of the 5th band near T Q and the 3rd and 4th bands 
near X is clearly correlated with, respectively, the A^^.j and A d xy characters. The 
cause for this lowering is therefore likely to be hybridization with the above-lying A 
bands. We saw in section 2 that the GdFeC>3-type distortion is larger in YTiC>3 than in 
LaTi0 3 , because Y 3+ is smaller than La 3+ , and because the Y 4c? bands are lower and 
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more narrow than the La 5d bands (figure 1). As a result, the A d character left behind 
to hybridize with the Ti t 2g band is different in the two cases, and this together with 
the electrostatics is the reason why the shape of these bands, and, hence, the orbital 
orders, to be considered below, are very different. 

In the far right-hand side of figure 12 we have shown the bandstructure calculated 
for the structure measured at 16 GPa, in which the JT distortion is nearly absent, the 
GdFe0 3 -type distortion is slightly increased, and the A positions have moved most [9]. 
As seen from table 7 and as mentioned at the end of section 3.5, at 16 GPa the crystal- 
field splitting t y z ° yz — t°® zxz 1S reduced from —230 to —167 meV, the crystal-field coupling 
^yz°xy~^xz°xy * s increased from —86 (see table 9) to —106 meV, and the anomalously small 
hopping integral, t^ yz , is decreased from —65 to —41 meV. As expected, the pseudo- 
gap and the A d xy character have both increased. So at 16 GPa, the bandstructure of 
YTiOs is distorted more, rather than less away from the cubic bandstructure of SrV03. 

Still, one may ask: is it the displacement or the scattering properties of the A 
ion by which the band shape (the pseudo-gap) is determined. To answer this, we 
have performed LMTO calculations of the bandstructure of LaTi0 3 using the crystal 
structure of YTi0 3 , both with the cell volume of YTi0 3 and with that of LaTi0 3 . We 
have also calculated the bandstructure of YTi0 3 using the crystal structure of LaTi0 3 , 
both with the cell volume of LaTi0 3 and with that of YTi0 3 . In all cases, the result 
was that, given the structure, it is the GdFe0 3 -distortion rather than the nature of the 
A ion which determines the shape of the t 2g band. 

Let us finally demonstrate that the JT distortion is not the reason for the 
development of the pseudo-gap. The blue t 2g bands in figure 13 are those of YTi0 3 
with the real, JT-distorted structure, which amounts to a 3% stretch of the 02 square 
into a rhomb in the y direction in subcells 1 and 3, and in the x direction in subcells 2 
and 4. The dashed red bands in the left-hand figure were obtained for the hypothetical 
structure with perfect Ti06 octahedra discussed in connection with table 7. We see 
that the band shapes, in particular those of the pseudo-gap, are the same, compared 
for instance with the huge difference in band shapes between LaTi0 3 and YTi0 3 . The 
results for the 16 GPa structure, in which the JT distortion is minimal, confirm this 
conclusion. 

3. 7. Crystal-field representation and orbital order 

Since for the titanates the calculated crystal-field splittings are an order of magnitude 
larger than the spin-orbit splitting (20 meV) and the magnetic ordering temperatures, 
it will prove useful to transform from the yz, zx, and xy Wannier orbitals to those linear 
combinations, |1) , |2) , and |3) , which diagonalize the on-site LDA t 2g Hamiltonian, t 000 , 
given in the first column of tables 1-5 and 7. Of these so-called crystal-field orbitals, |1) 
has the lowest energy and |3) the highest. In CaV0 3 the lowest orbital remains almost 
purely xy, in LaTi0 3 it is a fairly equal mixture of all three orbitals, and in YTi0 3 it is 
a mixture of yz and xy only. 
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Figure 12. LDA LMTO bandstructures of LaTi0 3 [12] and YTi0 3 at normal 
and 16 GPa pressure [9]. The bands have been decorated with A d xy and A 
(d xy + dy X + iz) /V3 partial-wave characters, averaged over all A sites. These 
characters are seen in figure 15 to be characteristic for the lowest crystal-field orbital. 



Specifically for CaVC>3, the eigenvalues relatively to Ep and the eigenvectors at site 
000 are respectively: 



ei £2 £3 
|1> |2> |3> 



\yz) \xz) \xy) 



For LaTi0 3 with the older data [44], we get: 
ei £2 e 3 



|1> |2> |3> 



\yz) \xz) \xy) 



meV 


and 




/ .202 


-.340 


.918 N 


.042 


.940 


.339 


v .979 


.030 


-.204 


meV 


and 




/ .715 


-.609 


.344 \ 


.353 


.739 


.574 


v .604 


.289 


-.743 i 



(17) 
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Figure 13. Left: YTi03 t% g bands (blue) for the real structure [20] and those obtained 
from the hypothetical structure without JT-distortion used also for producing table 7 
(dashed red). Right: YTi03 ti g bands for the real structure at normal pressure [9] 
(blue) and those calculated for the structure at 16 GPa [9], where the JT-distortion 
has essentially disappeared. The energy axis of the high-pressure bandstructure has 
been scaled such as to have the bands at normal and high pressure match well. The 
true scale may be found in figure 12. 



whereas with the recent LaTiC>3 data [12], results are somewhat different: 



ei e 2 63 



|1> |2> |3> 



\yz) \xz) \xy) 



meV 








/ .735 


-.538 


.413 
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.320 


.812 


.489 




v .599 


.227 


-.768 


/ 



(19) 



in particular for states |2) and |3) . For YTiC>3 the eigenvalues and eigenvectors are: 



ei e 2 e 3 



|1) |2> |3> 



289 488 620 meV 



\yz) \xz) \xy) 



( 


.781 


-.571 


.253 \ 




.073 


.319 


.945 


\ 


.620 


.757 


-.207 j 



(20) 



These crystal-field splittings are much larger than the spin-orbit splitting (~20 meV) 
and kT, and they agree with what was deduced (0.12-0.30 eV) from spin-polarized x-ray 
scattering for LaTiOs [47]. 

All three eigenf unctions and the two level splittings are shown in figure 14 for the 
titanates. The orbital of lowest energy, |1) , is of course the one for which the electrostatic 
attraction and the bonding hybridization with the A ions are maximized at the same 
time as the electrostatic repulsion and the antibonding hybridization with the oxygens 
is minimized. In LaTiOs, this orbital is directed towards the two La cations closest to 
Ti, the ones along the shortest diagonal. Also in YTiC>3 do the nearest cations attract 
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Figure 14. Ti t 2g crystal-field splittings and NMTO orbitals |1) , |2) , and |3) at site 

000 (subcell 1) for LaTi0 3 [12] and YTi0 3 [20] obtained by diagonalization of the 
on-site LDA Hamiltonian in the orthonormalized t2g-NMTO basis. With the older 
LaTi03 data [44], we get crystal-field levels of 0, 140, and 210 meV, whereas using 
the recent [12] 8K data yields 0, 207, and 221 meV. For the hypothetical structure of 
YTi03 without JT-distortion used in table 7, the crystal field levels are 0, 136, and 
259 meV. Orbital |1) has 97 (94) % of its charge density inside the central LagTiOe 
(YgTiOg) unit. The 3% difference is due to e s -character on the Ti 010 neighbours in 
YTiC>3. The same contour was chosen as in figure 8. The strong e g character on Ti 

001 for orbital |2) in YTiC>3 has been cropped in this figure. 

the red lobes of the lowest orbital, and for the same reasons, but in YTi03 there is also 
competition from the 2nd- nearest cations and, as a result, the blue lobes bond along the 
[111] diagonal with the Y 4d xy orbitals. The details of the hybridizations were explained 
for the constituent yz, xy, and xz orbitals in the preceding subsection. 

The difference between the lowest orbital in LaTi03 and YTiC>3 is brought out 
clearly in figure 15 where we have chosen a lower contour than in figure 14, as well as 
the orbital at site 010 (see also figure 16). In both materials does each red lobe bond 
to the closest-cation d 3z ^ _ x orbital directed along the shortest A-B-A diagonal, and 
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in both materials does each blue lobe form a complex 5-centre, B-2x02- Ad xy -01 
bond around the 2nd-nearest A-ion towards which the octahedron tilts. But the former 
mechanism prevails in La and the latter in Y titanate. The effect of A d 3z 2 i and A 
d xy hybridizations on the t<i g bandstructures were shown in figure 12. Finally, for the 
YTi03 orbital, we see that two oxygen p lobes which point opposite to the direction 
of their GdFe0 2 -type displacement, attain bonding e g character on their Ti neighbour. 
As was mentioned in subsection 3.2, this enhances the p lobe and bends it outwards. 
This hybridization -allowed only by the strong tilt- contributes as well to lowering the 
energy of orbital |1) . 

These quantitative differences in the shape of the lowest orbital in La and Y 
titanate are not caused by the small JT distortions. As an example, for YTi0 3 without 
JT distortion the on-site Hamiltonian in table 7 yields the following eigenvalues and 
eigenvectors: 



( e 1 e 2 e 3 ) = ( 341 477 600 ) meV 
( |1> |2> |3) ) = ( \yz) \xz) \xy) 



/ 


.766 


-.473 


.435 


\ 




.023 


.696 


.718 




\ 


.642 


.540 


-.544 


/ 



(21) 



For the lowest orbital, this is nearly identical with the results for the JT-distorted 
structure in equation (20), but it differs substantially for the higher orbitals, and that 
will prove important for the stability of ferromagnetic order in YTi0 3 (section 6). For 
the structure measured at 16 GPa, table 7 yields 

( ei e 2 e 3 ) = ( 365 584 726 ) meV 

/ .829 -.221 .514 \ ^ 
( |1) |2) |3) ) = ( \yz) \xz) \xy) ) -.243 .685 .689 . 

y .504 .693 -.514 y 

Here again, the shape of the lowest orbital is nearly the same as at normal pressure, 
although less cubic and more YTi0 3 'ish. This proves that the shape of the lowest 
orbital is determined mainly by the GdFe03-type rather than by the JT distortion. 

The shape of the orbital becomes crucial when it is placed in the crystal. As seen 
in figure 16, the orbital order is very different in the two titanates. In LaTi0 3 , the 
lowest orbital roughly has the be plane as mirror (x <-> y) , and this means that the 
glide-mirror operation from site 1 to 2 is roughly a translation, and so are therefore 
all cubic translations in the xy plane. In YTi0 3 , the lowest orbital does not have this 
symmetry, and the orbitals in the xy plane therefore avoid each other even more than 
in LaTi0 3 . As a result, the dominating integral, t\i° = t1\°, for hopping between two 
lowest orbitals on nearest neighbours parallel to the mirror plane is only half as large 
in YTi0 3 (—48 meV) as in LaTi0 3 (—98 meV). In both cases, this hopping integral 
is mainly via d xy character on the two A ions which are lst-nearest neighbours to one 
Ti and 2nd-nearest neighbours to the other. But whereas both large hopping integrals, 
tyl° yz and t° X y Xy , contribute in LaTi0 3 , only the latter does so in YTi0 3 . This may be 




Figure 15. Crystal-field orbital with the lowest energy, |1), at site 010 (subcell 2). 
A contour value 2/3 of the one used in the previous figures was chosen in order to 
exhibit the hybridization with the O p and A d 3z 2 and d xy states. In LaTi03, 
the hybridization with d^ z a _i on the two lst-nearest A ions dominates over the 
hybridization with d xy on the two 2nd-nearest A-ions towards which the octahedron 
tilts. In YTiC>3, the opposite is true. For the heavily distorted YTi03, there is also 
bonding from 02 p to e g on two Ti-neighbours. The A d Sz 2 _j and d xy hybridizations 
of the t2 g bands were shown in figure 12. 



seen in detail from the matrix transformation (23) and (24) given below. Both for La 
and Y titanate, there is only little 02-mediated hopping because the p character of one 
orbital hardly couples to the other orbital. 

The lst-nearest-neighbour hopping in the z direction, i^? 1 ; is antibonding (105 meV) 
in LaTiOs, with the sign-convention dictated by the mirror plane, and it is mainly via 
oxygen p. In fact, since according to the transformation (19), orbital |1) in LaTiC>3 is 
mainly composed of \yz) and \xy) , and since the latter hardly couples in the z direction, 
t?? 1 ~ — OM^yz and t yz x yz = —193 meV, as was discussed in connection with equation 
(5). In YTi0 3 , orbital |1) is almost exclusively composed of \yz) and \xy) , so that 
here we should have: t™ 1 ~ — 0.6t y ox yz . With 1^ yz anomalously small (—65 meV) due 
to coupling via Y Ad xy , t^? 1 becomes not only anomalously small, but even bonding 
(—38 meV) because now the small coupling (±64 meV) between \yz) and \xy) in the 
z direction cannot be neglected. This can be followed explicitly in equations (25) and 
(26) given below. 

Hence for the lowest crystal-field orbital the paths for hopping between nearest 
Ti neighbours in the xy plane are very different from those for hopping along z. Since 
exchange couplings are proportional to hopping integrals squared in conventional super- 
exchange theory, and since the spin-wave spectra in antiferromagnetic La as well as in 
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ferromagnetic Y titanate were both measured to be isotropic, this could be taken as 
an argument against the applicability of conventional theory for the low-temperature 
properties of these systems [40]. We shall return to this in section 6. 

Hopping integrals squared are very sensitive quantities, and since the residual 
cation covalencies in the materials with large GdFe0 3 -type distortions have completely 
ruined the cubic symmetry of the LDA t 2g Hamiltonians, let us explicitly perform 
the transformation of the dominating t w hopping integrals, t 010 , t 100 , and t 001 , from 
the cubic to the crystal-field representation in LaTi0 3 and YTi0 3 . We start from the 
cubic representation in tables 4 and 5, and use the transformation valid for site 000 
given by equations (19) and (20). When transforming the hopping integrals, we should 
remember that in the cubic representation (middle row in figure 8) those t 2g orbitals 
which would couple strongly if the crystal were cubic, carry the same name in all four 
subcells, whereas in the crystal-field representation, the orbitals are named and have 
signs following the space-group symmetry, that is, the bottom row in figure 8 and in 
figure 16. 

The t 010 matrix for hopping between crystal-field orbitals |1) , |2) , and |3) at sites 
000 and 010, which are related by the glide mirror (x <-> y) , is for LaTi0 3 : 



/ 



.735 
-.538 
.413 



.320 
.812 
.489 



.599 \ 
.227 
-.768 j 



( -185 
-43 



/ -98 -192 
4 76 
y 120 -40 



75 
1 

-33 

12 \ 
-22 
128 i 



-34 \ 

73 
-183 



/ 



/ .320 
.735 
v .599 



.812 
-.538 
.227 




meV 



meV, 



(23) 



with the order of rows and columns being like in equations (17)-(22), i.e.: yz, xz, xy 
for the cubic and 1, 2, 3 for the crystal-field representation. In the cubic representation, 
the elements along the diagonal of the hopping matrix are the well-known t n , t«s, t n . It 
may be realized that the hopping integrals, t^^,, resulting from this transformation are 
sums of contributions with varying signs. For YTi0 3 : 



/ 



V 



.781 
-.571 
.253 



-54 \ 
65 
-162 



/ _ 



V 



073 


.319 


.945 
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781 


-.571 


.253 




620 


.757 


-.207 


) 



meV 



.073 .620 \ / -184 70 
.319 .757 -41 28 

.945 -.207 J \ -64 -22 

/ -48 -191 -130 \ 
-84 -13 44 
94 11 -73 

Note that at the second site 010, we must exchange the yz and xz orbitals. For the t 100 
matrix of hopping integrals between the orbitals at sites 000 and site 100 = 010 + a the 
transformation is the same, and we simply get: 



V 



meV. 



(24) 



/ 



-100 _ j-U 
mm' m'm' 



-010 



For these hoppings parallel to the mirror plane, the differences between La and Y- 
titanate are rooted more in the different crystal-field eigenvectors than in the different 
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hoppings between the yz, xz, and xy orbitals. 

The t 001 matrix for hopping between crystal-field orbitals |1) , |2) , and |3) at sites 
000 and 001, which are related by the A-Ol mirror plane perpendicular to the z axis 

/ -193 -42 -59 \ / 
-42 -208 52 
59 -52 -22 



—z), is for LaTi0 3 : 
.735 .320 .599 \ 
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meV, 



(25) 



and for YTi0 3 it is: 

/ .781 -.073 .620 \ / 
-.571 .319 .757 
.253 .945 -.207 

/ -38 -21 
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meV 



(26) 



Here we must flip the sign of the yz and xz orbitals at site 001. In case of hopping 
perpendicular to the mirror plane, not only the different crystal-field eigenvectors 
matter, but also the large non-cubic perturbations of the hopping integrals in YTiC>3. 

In the cubic representation, the matrices of lst-nearest-neighbour hopping integrals 
have two large diagonal elements, t n , and all other (inter-orbital) elements are small, 
except in YTi0 3 , where tfj£ yz = —65 meV is anomalously small and the inter-orbital 
elements are of similar size. In the crystal-field representation, the matrices of lst- 
nearest-neighbour hopping integrals do not have this form at all. In particular in YTiOs, 
the hopping between orbitals |1) is anomalously small and, except in one case, smaller 
than the hopping from orbital |1) to orbitals |2) and |3) . As we shall see in section 6, 
this is what makes YTi03 ferromagnetic at low temparature. 

In the crystal-field representation, we do not mirror the orbitals like in equation (4) 
in order to obtain the values of integrals corresponding to hops starting from another 
site than the first. Here the rules are simply: 



t 



(0,0,l),(s,j/,*+l) 
m',m 



,xyz 
m'm 
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(0,1,1), (x,y+l,z+l) 
m' ,m 
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A0,0,0),(x,y,z) _ ,xyz 
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In figure 17 we show the bandstructures decorated with B t2 g partial- wave characters 
in the cubic and the crystal-field representations at site 000 or, equivalently, according 
to the bottom row of figure 8 (see also figure 16). In general, the mixing of characters 
is considerable in both representations, although the crystal-field splittings - which are 
small on the scale of the bandwidth - do cause a slight preference for lower-lying bands 
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Figure 16. Orbital of the lowest energy, |1) , placed in the 4 subcells (1=000, 2=010, 
3=001, and 4=011) according to the mirror and glide-mirror symmetries, i.e. like in 
the bottom row of figure 8. This illustrates the orbital order which will crystalize out 
due to the Coulomb correlations (see subsection 5.3). 

to have stronger |1) character (red) and higher-lying bands to have stronger |3) character 
(blue). Nevertheless, there are two important cases of character separation: the bands 
at the lower edge of the pseudo-gap, which move down and get increasingly occupied as 
we proceed along the series, have strong |1) character, and the 4th band near r o , which 
moves upwards and almost empties, has strong |3) character. As we have seen in figures 
12 and 15, the |1) bands are those with strong A d 3z 2^_ 1 and A d xy characters. The 
dominant role of those bands in the development of the pseudo-gap is also evident from 
the top right and the bottom row parts of figure 18, exhibiting the on-site 000 elements, 

N mm i (e) = Mr "M ( k ) 5[e-Ei (k)] u* nm , 4 (k) , 

of the density-of-states (DOS) matrix. Here, Uj (k) is an eigenvector of the LDA 
Hamiltonian (7) in a t^g Wannier-function basis. 

4. Multi-band Hubbard Hamiltonian and its solution in the dynamical 
mean-field approximation 



Having discussed at length the one-electron, "chemical" part of the Hamiltonian, 
H^ DA , we now add the on-site Coulomb repulsion and obtain the multi-band Hubbard 
Hamiltonian: 



H t 2g 



n 



LDA 

*2o 



R 



(27) 




Figure 17. LDA LMTO bandstructures of the orthorhombic perovskites decorated 
with B partial-wave characters at site 000 in the cubic and the crystal-field 
representations, or equivalently, according to the bottom row of figure 8 and figure 16. 
For CaVC>3 |1) « \xy) , \2) \xz) , and |3) » \yz) . The points of high symmetry are 
those mentioned explicitly in the caption to figure 9, and not any of their equivalents. 
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Figure 18. The LDA on-site DOS matrices N mm i (e) for site 000 (subcell 1) in the 
representation of the ti g Wannier functions. The zero of energy is the Fermi level. The 
DOS unit is states/(eV-spin-AB03). The insets show the off-diagonal elements on a 
reduced energy scale, but the same DOS scale. Top two rows: cubic basis as defined 
in the bottom row of figure 8. Bottom row: crystal-field basis as defined in figure 14 
and equations (17)-(20). The tick marks indicate the positions of the t^ g crystal-field 
levels. The occupancy of the lowest crystal-field orbital is 0.45 in LaTiOa and 0.50 in 
YTi03. The Fermi levels for occupations with one and two electrons are marked; the 
latter is the majority-spin Fermi level for a d 1 ferromagnet in the Stoner model. The 
older data [44] was used for LaTi03. 
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Here, U is such that for a doubly occupied site {fig) > its 15 degenerate states are as 
follows: a 9-fold degenerate triplet with energy U — 3 J, a 5-fold degenerate singlet with 
energy U — J, and a singlet with energy U + 2 J. For the analogous ^-configuration these 
states would be 3 P, 1 D, and 1 S, respectively The average d 2 energy is 

U - 2 J = U'. (28) 

U is the Coulomb repulsion between two electrons in the same orbital, and J is the 
Hund's-rule exchange coupling. 

Since TCt 2g involves only the correlated t 2g orbitals, the number of correlated 
electrons is fixed and the double-counting correction therefore amounts to an irrelevant 
shift of the chemical potential, which we omit in (27). For transitions from or to a 
non-correlated LDA band, such as those seen in figures 1, 3, and 4, there is a double- 
counting correction, but it is presumably < 1 eV, because the number of B t 2g electrons 
in the LDA calculation is not far from 1. In those figures, the LDA t 2g band should 
therefore simply be substituted by the correlated t 2g band with the Fermi levels lined 
up. 

In general, an LDA+?7 Hamiltonian like (27) depends on which electronic degrees 
of freedom, i.e. orbitals, are included, and how the correlated orbitals are chosen. This 
means that in order to be able to treat electronic correlations beyond the LDA, we 
have to depart from the ah initio philosophy and make a system- and even property- 
dependent choice of an appropriate low-energy Hubbard Hamiltonian. In the present 
work for the d 1 perovskites, the correlated orbitals are taken as the set of localized t 2g 
LDA Wannier functions. The degrees of freedom from all other bands are neglected, not 
only in the correlation term, but also in H LDA . For instance, the B e g band, whose centre 
is 3 eV above that of the t 2g band (see figures 3 and 4), is not treated as correlated, but 
the e g characters are downfolded into the t 2g orbitals as seen in figures 14 and 15. For 
the Mott transition, this seems to be the appropriate treatment, as is also suggested 
by the work of Manini et al. [15]. For super-exchange couplings this also seems to 
hold, as we shall argue in section 6, but for high-energy optical or inverse-photoemission 
spectra, it is clearly inappropriate. Here, e g degrees of freedom must be taken explicitly 
into account in the Hubbard Hamiltonian. This can be done by using as basis either 
the present t 2g basis augmented by the Wannier functions for the LDA e g band, or the 
set of Wannier functions for the entire d (t 2g + e g ) band. The advantage of the former 
basis set is that it has no single-particle coupling between the t 2g and e g bands, and the 
advantage of the latter is that it is more localized. Anyhow, in the present work we do 
not compute high-energy properties. 

In principle, the on-site Coulomb matrix should be calculated from the Wannier- 
function basis used, but in the present case it was simpler to keep a single, adjustable 
parameter, U, and - as we shall see - this yields more insight. The value of this 
parameter we took to be the same for all four materials. For the Hund's-rule coupling 
we used the following: 

J = 0.68 eV and 0.64 eV, 
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for the vanadates and titanates, respectively. These are atomic Hartree-Fock values 
times 0.8 in order to account for the screening in the solid [25]. 

The Hubbard Hamiltonian (27) is solved in the DMFT [35], i.e. under the 
assumption that the elements, £_R m ,R' m ' (u) , of the self-energy matrix between different 
sites can be neglected. The self-energy is thus assumed to be an effective, energy- 
dependent and complex crystal-field term. In this case, the on-site Green matrix is: 

G Rm , Rm/ (u) = V { [u - E (u) - H LDA (k)] _1 | , (29) 

— I J Rm,Rm' 

k 

and the many-body lattice problem is then mapped onto an Anderson impurity problem 
in which the inverse of the bath Green function of the uncorrelated host is required to be 
G (uj)^ 1 + £ (u) . Solution of the Anderson impurity problem must now yield the same 
GRm, Rm ' {w) as equation (29), and this is a self-consistency condition for determination 
of the self-energy matrix, S^ mj ^ m / (u) , and the host. 

The bottleneck in a DMFT calculation is to solve the Anderson impurity problem, 
and only a few correlated orbitals can be handled at the moment. In the orthorhombic 
perovskites all 4 B-sites are equivalent, so the t 2g impurity problem involves only 3 
correlated orbitals. To solve it, we employed the numerically exact Hirsch-Fye [48] 
quantum Monte Carlo (QMC) method. In order to access temperatures down to 770 K, 
we used up to 100 slices in imaginary time and about 10 6 QMC sweeps. To reach 
convergence, 15-20 DMFT iterations were needed. Finally, the spectral function was 
obtained on the real cj-axis by analytical continuation using the maximum entropy 
method [49]. Unfortunately, this does not provide us with the self-energy matrix for 
real u, so we do not obtain correlated bandstructures in the present applications. 

In order to be able to perform the DMFT calculations we adopted the usual 
approximation of keeping only the density- density terms in the on-site Coulomb- 
repulsion, i.e. 

{U if m — m! and a ^ a' 

U - 2 J if m^m' and a ^ a' 
U — 3J if m 7^ m' and a = a' 

Here, nR ma = o} Rmcr aR mcF , and a} Rma creates an electron with spin a in a localized orbital 
m at B-site R. In the summation, at least two of the indices on the operators must be 
different. The repulsion averaged over all doubly occupied states remains as in equation 

(28) [50]. 

In previous implementations of the LDA+DMFT method [36, 37] it was assumed 
that the on-site block of the single-particle Green function is diagonal with identical 
elements in the space of the correlated orbitals, and the latter were usually taken as 
orthonormal LMTOs [41], approximated by truncated and renormalized partial waves. 
These approximations mean that the partial waves not belonging to the irreducible 
representation of the correlated LMTOs -and therefore arising from the tails of LMTOs 
on neighbouring sites- are neglected. As a consequence, the self-energy in equation 

(29) merely distorts the energy scale of the LDA DOS, from which the on-site G (e) 
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can then be obtained by Hilbert transformation. Although this may be appropriate for 
e.g. cubic systems, it is clearly not appropriate for the series of materials considered 
in the present paper. Here, the self-energy must be allowed to be an effective crystal 
field - enhanced, energy dependent, and complex through the influence of Ur - which 
produces k-dependent broadenings and distortions of the bandstructure of the kind we 
have been discussing extensively in section 3.6. Specifically, in table 9 we should perform 
the substitutions: 

^000 J)00 , y /A /on\ 

L XiXj,XjX k L XiXj,XjX k ' ^XiXj,XjX k \ oyj J 

It is not correct to merely distort the energy scales of the different partial DOS functions 
differently. 

In our new implementation [21] of the LDA+DMFT method, the highly accurate 
NMTO method is used to create a localized set of correlated orbitals, which then defines 
the Hubbard Hamiltonian. We choose the localized orbitals to span a subset of LDA 
bands around the Fermi level exactly -the t 2g bands in the present application,- and 
we choose to orthonormalize them symmetrically so that they become a set of Wannier 
functions. With this set, the matrix elements are then calculated correctly (except the 
matrix elements of the on-site Coulomb repulsion). Hence, we take all components of 
the self-energy matrix S mm / (u) between different orbitals on a given site into account. 
From this 3x3 matrix we use the Pbnm symmetry (figure 2) to construct a 12x12 
block-diagonal self-energy matrix, which is then inserted in expression (29). The 12x12 
matrix uj — £ (cu) — H LDA (k) is now inverted as a function of k, and, finally, the on-site 
3x3 block is summed over k to yield the 3x3 on-site G (e) matrix. 

A recent LDA+DMFT calculation for Lai^Sr^TiOa [39] used the crystal-field 
representation and neglected the off-diagonal elements of the DOS and the self-energy. 
As seen in the bottom row of figure 18, this may be an intelligible approximation. 
Nevertheless, even with the crystal-field eigenvectors frozen at the LDA values, the 
S (u) + H LD (k) bandstructure should be allowed to change as the effective crystal- 
field splittings change from iteration to iteration. 

Another LDA+DMFT implementation based on Wannier functions was proposed 
after the completion of the present work [51]. Applications to cubic SrVOs using t 2g 
Wannier functions yielded low-energy spectral functions similar to those presented in 
the next section. It is not obvious to us how t 2g Wannier functions can provide matrix 
elements of the self-energy connecting to other bands, except for the trivial double- 
counting shifts mentioned above. A good idea in that paper is to use the on-site Green 
function (29) to evaluate the contribution to the charge density from the correlated 
orbitals and, hence, to adjust the Wannier functions self-consistently, away from their 
LDA values. For the t\ g titanates, a useful application would be to calculate the JT 
distortions from the LDA+DMFT charge density. 
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Figure 19. DMFT spectral function for [7=3.5 (red), 4.0 (bold red), 4.5 (green), 4.75 
(bold green), 5.0 (blue), and 6.0 (bold blue) eV. The average Coulomb repulsion is 
U'=U - 1.28 eV. kT= 100meV-1200K. For LaTi0 3 , the older structural data were 
used [44]. Extensive DMFT calculations were not performed with the recent data [12]. 

5. High-temperature properties 

5.1. The Mott transition 

With LDA+DMFT we computed the spectral functions for all four materials as a 
function of U and at kT=100 meV~1200K, i.e. way above any magnetic ordering 
temperature. We do not investigate the complete phase diagram or the order of the 
Mott transition; rather, we estimate the critical value, U c , for the high-temperature 
region where the system is paramagnetic and for which we find a unique solution of the 
DMFT equations. 

For the titanates, the spectra are shown in figure 19. We found that the materials 
become insulating when U' > U' c with the values of U' c given in the first row of table 
10. Cubic SrVOa remained metallic at the highest value of U' considered (4.6 eV), and 
the calculation employed the older structural data [44] for LaTiO"3. The decrease of the 
critical value seen along the series is dramatic, and if we form the ratio with the rms 
bandwidth W from table 8, we find the numbers given in the central row of table 10. 
This ratio is not constant, but decreases along the series by 35% from CaVC>3 to YTiO"3. 
Therefore, the Mott transition is not driven merely by the narrowing of the ti g band, 
but also by its deformation. That is unexpected. 

The importance of orbital degeneracy for the Mott transition was first pointed 
out by Gunnarsson et al. [14], who argued that due to the increase in the number 
of hopping processes in many-body theory compared with band theory, there is an 
enhancement of the hopping, which for a half-full band is approximately proportional 
to the square root of its degeneracy. The value, U' c , necessary to cause a Mott 
transition for a given bare bandwidth, W, therefore increases with the degeneracy. 
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Table 10. Critical Coulomb repulsion, U' c , in eV and relatively to the rms bandwidth, 
W. The crystal- field splitting between the lowest levels is Ai 2 . 





SrV0 3 [42] 


CaV0 3 [43] 


LaTi0 3 [44] 


LaTi0 3 ([12]) 


YTiO 3 [20] 




> 4.6 


4.4 


3.6 




2.5 


U'JW 


(1.96) 


1.84 


1.65 




1.34 


A 12 /W 





0.031 


0.064 


0.093 


0.106 



Imagine now with Manini et al. [15] that this degeneracy is split by a small crystal 
field: upon increasing U', and thereby reducing the width of the quasiparticle peak 
towards zero, ZW ~ a [1 - U'/U' c (N)] W, for some value of U' < U' c (N) the reduced 
bandwidth will have reached the same size as the crystal-field splitting. At that 
point the degeneracy is effectively decreased from N to N — n, and herewith U' c , so 
that a Mott transition may be triggered. One therefore expects a critical crystal-field 
splitting, A c ~ a [1 — U' c (N — n) jU' c (JV)] W. Manini et al. specifically solved a simple 
Hubbard model with one electron in two identical, non-interacting bands whose on- 
site energies differed by A For A = 0, they found U' c (2) jW ~ 1.8, and for A — > oo, 
they found U' c (l) /W ~ 1.35, so that U'J2) /U' c (l) ~ 1.33. The surprising result was 
that U'JW decreases from the —filled two-band limit to the —filled one-band limit 
for A increasing from to merely 10% of the bandwidth, i.e. A c /W ~ 0.1, so that 
a ~ 0.1/ (1 — 1.35/1.8) ~ 0.4. The crystal field splitting does, however, only control 
U'JW for the Mott transition within the limits from U' c (2) /W to U' c (l) /W; if for 
A > A c , U' /W is increased from below 1.35, then the system undergoes a phase 
transition from a two-band metal to a one-band metal before it goes from a one-band 
metal to a one-band Mott insulator at U'JW = 1.35. 

For our case of a |-filled three-band system, Koch et al. [52] solved a Hubbard model 
for the ti u -band in Cgo and found that hopping is enhanced by merely ~ (l + v^2) /2 = 
1.21. Thus, the degeneracy reduction expected from the crystal-field splitting is at most 
30%. The striking fact about our values of U'JW is that they decrease by considerably 
more: if we use Manini et a/.'s a- value, then the 3% crystal-field splitting in CaV0 3 
has reduced U'JW by about 6% from U' c (3) /W, the value appropriate for SrV0 3 and 
listed in parentheses in table 10. Our LDA+DMFT calculations thus yield a decrease 
in U'JW of ~50% when going from cubic SrV0 3 to YTi0 3 . That is not possible unless 
the width of the lowest subband narrows more than the rms bandwidth, that is, unless 
the band deforms along the series, an effect not considered previously. We shall return 
to this point in section 7. 




Figure 20. DMFT spectral functions for U=5 eV and T=770 K (bold) and LDA 
total density of states (weak). For LaTi03, the old structure was used. Reproduced 
from [21]. 



5.2. Spectral functions for U —5 eV and comparison with experiments 

The main features of the photoemission spectra for all four materials, as well as the 
correct values of the Mott-Hubbard gap for the insulators, are reproduced by taking 
U material independent and equal to 5 eV. This is satisfying because U is expected to 
be similar for vanadates and titanates, or maybe slightly smaller for the latter [25, 53]: 
Bocquet et al. fitted 2p core-level photoemission spectra with a model Hamiltonian 
containing all Op and Bd orbitals and obtained U' xs 4.0 eV for La and Y titanate [54]. 
This is consistent with the value U' — 5.0 — 2 x 0.64 = 3.7 eV used by us for the less 
localized t 2g orbitals. Similar values of U' have been used by other authors [38, 51] 

In Figure 20 we show the spectral functions together with the total LDA DOS. The 
vanadates exhibit a quasiparticle peak and are therefore metallic, while the titanates 
are Mott insulators. 

For cubic SrV03 we reproduce the results of previous calculations [38, 55]: the 
lower Hubbard band is around — 1.8 eV and the upper Hubbard band around 3eV. For 
CaV0 3 , the quasiparticle peak loses weight to the lower Hubbard band, which remains 
at — 1.8 eV, while the upper band moves down to 2.5 eV. These results are in very good 
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agreement with photoemission data [38, 56]. From the linear regime of the self-energy 
at small Matsubara frequencies we estimate the quasi-particle weight to be Z = 0.45 
for SrVC>3 and 0.29 for CaVC>3. For a k independent self-energy, as assumed in DMFT, 
this yields m*/m = 1/Z = 2.2 for SrVC>3 and 3.5 for CaVC>3. These results are in good 
agreement with the optical-conductivity values of 2.7 and 3.6 [4], as well as with the 
values 2-3 obtained by thermodynamics [5, 6] and dHvA experiments [46]; they show 
that SrV0 3 and CaV0 3 are rather similar, with the latter slightly more correlated than 
the former. Similar conclusions were drawn in [38]. 

For the titanates, the lower Hubbard band is around — 1.5 eV, in accord with 
photoemission [54, 53]. Despite very similar bandwidths, the gaps are very different, 
approximately 0.3 eV for LaTi0 3 and approximately leV for YTi0 3 , and this also 
appears to agree with optical experiments [7]. 

5.3. Orbital polarization 

Diagonalization of the matrix of occupation numbers, 

1 f 

obtained with LDA+DMFT for values of U' exceeding the t<i g bandwidth, W, reveals 
that - within the numerical accuracy - the eigenvectors, w XiXj>m , are the same as those 
which diagonalize the on-site LDA Hamiltonian matrix, that is, n and t° commute. For 
the titanates, the eigenvectors are the ones given by equations (18), (19), and (20), 
and the corresponding eigenfunctions are the crystal-field orbitals shown in subcell 1 
in figure 14. However, upon increasing U the orbital polarizations increase around the 
metal-insulator transition -from the small LDA values given in the top and bottom 
rows of figure 18- and become nearly complete after the Mott transition. For this to be 
true, in the insulators the DMFT must create a self-energy matrix which strongly and 
uniformly enhances the non-cubic part of the on-site LDA Hamiltonian matrix, i.e. for 
E in equation (30) we must have: 

^xiXj,xjx k ~ y ] w XiXjtm n m U m {uj) w XjXktm , (31) 

m 

so that in the crystal-field representation 

E mm , {uj) « 5 mm 'S m {uj) = o~ mm >n m U' m {uj) , (32) 

where, as U' increases beyond U' c , n m — > ~ 1 while Re U' m {uj) tends towards a function 
which is roughly - — |[/' + 0{uj) for uj negative and - + \U' + 0(u) for uj positive, 
and Re U' 2 {uj) and Re {uj) tend towards functions which are roughly - + |[/' + 0{uj) 
independently of the sign of uj. Here, a is a constant. This would mean that the Mott 
transition takes place in only one band, and that this band is the one described by the 
crystal-field orbital with lowest energy. Static mean-field methods like Hartree-Fock and 
LDA+f/ have no ^-dependence of U' m {uj) , and are therefore unable to split the half-full, 
lowest band into two. Unfortunately, we cannot present figures of E m {uj) because the 
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Figure 21. Diagonal elements of the spectral function in the crystal-field basis 
calculated with LDA+DMFT (17=5 and fcT=0.1 eV). The corresponding figure for 
C/=0 and kT=0 is the bottom part of figure 18, which also has the same colour 
coding: |1) red, |2) green, and |3) blue. The occupation of |1) is 0.88 for LaTi03 
(older structure) and 0.96 for YTiC>3. Below Ef there is one, and above there are five 
electrons. The off-diagonal elements are two orders of magnitude smaller than the 
diagonal ones. 



methods we use to solve the Anderson impurity problem do not enable us to evaluate 
the self-energy for real uj. 

Figure 21 shows the DOS for the titanates calculated with LDA+DMFT (£7=5 eV 
and kT=0.1 eV) and projected onto the common eigenvectors of the occupation numbers 
and the crystal-field Hamiltonian. Hence, for LaTi03 the lowest crystal-field orbital 
contains 0.49 [12] (0.45 [44]) electrons when £7=0 and 0.91 [12] (0.88 [44]) electrons when 
£7= 5 eV. For YTi03, the respective occupations are 0.50 and 0.96. The nearly complete 
orbital polarization found for the two insulators indicates that correlation effects in 
the paramagnetic Mott insulating state considerably decrease orbital fluctuations, and 
makes it unlikely that YTiC>3 is a realization of an orbital liquid [57] . In LaTi03 some 
orbital fluctuations are still active, although quite weak. 

Most recent NMR [58] and x-ray scattering [59] experiments point to an orbital 
order which is very similar to the one obtained in our calculations. For LaTiOs our 
occupied orbital is very similar to the one arrived at slightly earlier by Mochizuki and 
Imada [24] through model calculations to which we shall return in more detail in section 
6. Also some of the results in the recent LDA+DMFT calculation for Lai-^Sr^TiOs by 
Craco et al. [39] mentioned at the end of section 4 are similar to ours. 

Our occupied orbital for YTi03 is similar to the one obtained with the spin- 
unrestricted LDA (actually GGA) scheme for T=0 [32] . This is understandable because 
the LDA properly accounts for the covalency and the electrostatics, and because the 
occupancy in the spin-unrestricted GGA is very similar to the one obtained by occupying 
the LDA DOS in the bottom right of figure 18 with a single spin-up electron, i.e. filling 
it to the Fermi level corresponding to n=2. The width of the band filled with one 
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spin-up electron is seen to be 0.8 eV so that in order to stabilize this ferromagnetic 
state, a Hund's-rule coupling constant exceeding 0.4 eV is required; this is well satisfied 
by our J=0.64 eV. In being metallic, this ferromagnetic ground state is, however, 
incorrect and the LDA+[7 scheme was needed to produce an insulator [32]. Our 
orbital in equation (20) agrees almost perfectly with the one deduced from NMR [33]: 
|1) = 0.8 \yz) + 0.6 \xy) , and well with the one deduced from neutron scattering [34] and 
a small correction [60]: |1) = 0.7 \yz) + 0.7 \xy) . Also resonant a>ray scattering confirms 
the orbital order [60, 61]. Mochizuki and Imada obtained a similar orbital to which we 
shall return in more detail in section 6. 

To visualize the orbital order, the nearly full orbital was placed on each of the 4 Ti 
sites in figure 16. Despite the fact that LaTiOs and YTiOs have the same space group, 
the orbital orders look very different. This difference is, however, quantitative rather 
than qualitative; as mentioned in section 3.7, it reflects the extent to which the orbital 
has the be plane as mirror. We emphasize once more that the two different types of JT 
distortions observed in LaTiOs and YTiC>3 are not the cause for the difference in the 
orbital orders. This was discussed for YTi0 3 in section 3.5, was clearly shown in figures 
12 and 13, and was finally explicitly brought out in section 3.7 by the eigenvectors (21) 
and (22). It is however obvious that the oxygen octahedron will relax to the shape of the 
charge density of the localized electron, and this seems to be the reason why in YTiC>3 
the 02 square is stretched by 3% into a rhomb along y in subcells 1 and 3, and along 
x in subcells 2 and 4 [20], while in LaTiOs, it is stretched into a rectangle by 3% along 
a [12]. This is remarkable because in systems such as LaMnOs, where the low-energy 
bands are of e g type, orbital order and large (10%) JT distortions occur together, and 
the different types of orbital orderings are therefore often classified according to the type 
and spatial arrangement of the JT distortions. We have seen that this does not apply 
when the low-energy bands are of t 2g type, first of all because the t 2g orbitals have pdn 
coupling to the same oxygen orbitals as the A ions have pda coupling to and, secondly, 
because this pdn coupling of the t 2g electrons is much weaker than the pda coupling of 
the e g electrons (see figures 5 and 6). In conclusion, the crystal-field splitting in the t\ g 
perovskites is due to the GdFeOs-type distortion. 

We have thus seen that it is the Coulomb repulsion which causes the electron to 
localize, but it is chemistry which sets the stage: it determines U c and selects the orbital 
to be occupied. 

5.4- Pressure-induced metallization 

Loa et al. recently performed high-pressure experiments on the titanates, in which they 
monitored the structures and the optical gaps [9]. Some of their results are shown 
in the left-hand side of figure 22. At HGPa, LaTi0 3 becomes metallic, its volume 
contracts slightly, but no apparent change of internal parameters occurs. For YTi0 3 , 
the optical gap decreases with increasing pressure, but the material remains insulating 
up to at least 17GPa. By extrapolation of this pressure dependence, it was estimated 
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Angle(°) P(GPa) E Gap (eV) Bandwidth W(eV) 

Figure 22. High-pressure metallization of LaTi03 (blue) and YTiC>3 (red). 
Experimental [9] and theoretical data (horizontal axes) as a function of the reduction of 
volume from that at normal pressure (vertical axis, log scale). 9 and <p are the tilt and 
rotation angles. Rietveld refinement of the crystallographic data was only possible for 
V > 0.93Vb- For YTi03 the Mott transition has not been observed within the available 
pressure range, and the dashed line extrapolating the optical gap to zero is a theoretical 
guess from the right-hand side of the figure. For LaTiC>3, the Mott transition takes 
place at 11 GPa and y=0.94Vb, where a discontinuous contraction is observed. W 
is the calculated LDA ^g-bandwidth (distance between band edges), which is seen to 
follow a power law in the Ti-Ti distance d oc V 1 ^ 3 . Since the experimental GdFeC>3- 
type distortion seems to be independent of volume, at least up to V=0.9Vq, the values 
[7^=3.6 and 2.5 eV (table 10) obtained from our LDA+DMFT study of the Mott 
transition (figure 19) can be used to predict the metallization volume. 

that metallization will occur above 40 GPa [9]. 

We may compute the metallization pressures using LDA+DMFT, provided that we 
know the structure. Luckily, the experimental data in the left-hand side of figure 22 give 
no reason to expect that the structures change by anything but a uniform compression; 
the GdFeC>3-type distortion seems to be fairly constant both for LaTi03 and for YTiOs, 
at least up to a volume reduction of 0.9. We can therefore use our previous LDA+DMFT 
calculation of U' c , the value of the average on-site Coulomb repulsion necessary to drive 
the metal through a Mott transition (figure 19 and table 10), to predict the LDA t2 9 
bandwidth necessary to drive the insulator metallic. This gives: 

5 — 2 x 64 

W c ~ 2.09 eVx — = 2.2 eV for LaTi0 3 

3.6 
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and 

W c ~ 2.05 eVx 5 - 2 X 0,64 = 3.0 eV for YTi0 3 . 
2.5 

Since only the width -not the shape- of the band changes, we have used here the 
simplest measure of the bandwidth, namely the energy distance between the band edges 
(W t2g ). To find the critical volume reduction, V c /V , we now merely need to compute 
this bandwidth as a function of volume. The result is shown in the right-hand side of 
figure 22. The critical volume and corresponding pressure is then 

V c ~ 0.96Vb and P c ~ 7 GPa for LaTi0 3 , 

which is in reasonable accord with the experiments, and 

V c ~ 0.78 V and P c ~ 100 GPa for YTi0 3 . 

For YTiO"3, V c is uncertain because the GdFe03-distortion is only known to be constant 
for volumes larger than 0.9Vo, but the corresponding, sketched behaviour of the optical 
gap does not look unreasonable. The critical pressure is more uncertain than the volume 
because it hinges on an extrapolation of the pressure-volume relation as well. 

A most recent, refined study of YTi0 3 [9] revealed that between 9 and 14 GPa the 
JT distortion essentially disappears, but the distortion of the Y sublattice increases. 
This was mentioned at the end of section 3.5, and the 16 GPa structure was used to 
calculate the bandstructure shown on the right-hand sides of figures 12 and 13. This 
bandstructure justifies our approximation of assuming a rigid bandshape, although the Y 
positions shift so as to increase the residual covalency and, hence, to increase the pseudo- 
gap. A corresponding correction of U c would be towards a slightly lower value, and that 
would increase our estimate -to above 100 GPa- for the pressure where metallization 
should occur if no further structural changes were to take place. Most importantly, 
the fact that the JT distortion is strongly reduced in the 16 GPa structure does not 
mean that the orbital order has changed drastically. That was explained at the end 
of the preceding subsection. The reduced JT distortion will, however, change the low- 
temperature magnetic order from ferromagnetic to A-type antiferromagnetic, as well as 
the ordering temperature. This we shall show in section 6. 

5. 5. Onset of optical absorption in the titanates 

Let us calculate how the difference between LaTi0 3 and YTi0 3 shows up in the 
optical transitions between the lower and upper t 2g Hubbard bands seen in figure 21. 
Experimentally, the onset of the optical conductivity is more gradual in LaTi0 3 than 
in YTi0 3 [9, 62, 63]. 

The optical conductivity has previously been evaluated for LaTi0 3 in LDA+£7 [64] 
and LDA+DMFT [65] calculations. However, the GdFe0 3 -type distortion, crucial to 
the present study, was neglected in both calculations where, on the other hand, the 
optical matrix elements were treated with more care than in our study. 
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Figure 23. Real part of the optical conductivity obtained by LDA+DMFT with U=5 
eV, kT— 0.1 eV, and treating the matrix elements as explained in the text. The recent 
structural data for LaTiC>3 was used [12]. 



The optical conductivity can be expressed in terms of the current-current correlation 
function [66] . In the dynamical mean- field approximation there are no vertex corrections 
and the conductivity can therefore be written as [67, 68, 69]: 




Here, n(u) is the Fermi distribution function, My (k) the optical matrix element, and 
Aij (k, u) = -Im { \u - S (u) - H LDA (k)l _1 } 

7T L i id 

the spectral function at real energy. The latter can be obtained by analytical 
continuation from the QMC data created in the course of an LDA+DMFT calculation. 
We shall assume that all transitions between t 2g crystal-field orbitals (see figure 21), 
and allowed by the Fermi functions, have equal probability. Compared with charge- 
transfer transitions (0 2p-Ti3c? and 2p-Ad), these transitions are all weak, and we 
do not expect some to be stronger than others. With this approximation, equation (33) 
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reduces to: 

cr(u) oc — ^ / dJ [n{uj') — n(u + a/)] ^ A^(k, uj') ^ Ajj(k, uj' + uj), 

k ij ij 

where i and j now run over the crystal-field orbitals. 

The real part of the optical conductivity obtained from our LDA+DMFT (U =5 eV, 
kT=0.1 eV) calculations are given in figure 23. By going from the indirect DOS gaps 
in figures 20 and 21 to the direct gaps in figure 23, the latter are seen to be increased 
by approximately 0.5 eV, to 0.8 eV in LaTi03 and to 1.5 eV YTi0 3 . We observe that, 
even without detailed matrix-element effects, the onset of inter-band transitions is more 
gradual in LaTi0 3 than in YTi0 3 . That is presumably the same trend as seen for the 
upper Hubbard band in figure 21, where the lower edge of bands |2) and |3) rises more 
gradually in LaTi0 3 than in YTi0 3 ; in LaTi0 3 the peak is reached after 1.5 eV, but in 
YTi0 3 after 1.0 eV. 

6. Low-temperature, magnetic properties of the titanates 

6.1. Introduction 

The origin of the magnetic phases of LaTi0 3 and YTi0 3 have been a puzzle for decades 
[2, 3]. Whereas LaTi0 3 is a G-type antiferromagnet with T/v=150 K and a small moment 
of 0.57/is, YTi0 3 is a ferromagnet with a low Curie temperature, Tc=30 K, and a good- 
sized moment, 0.8/zb [11, 12]. 

The early idea was that in Mott insulators, the magnetic structure is rooted in 
orbital order. Hence, the different magnetic ground states in LaTi0 3 and YTi0 3 ought 
to arise from different orbital orders. Orbital order, presumably caused by electron 
correlation, should be accompanied by a JT elongation of the TiO@ octahedron, and 
in perovskites the stacking along the c axis was found to be either parallel (d type) 
or perpendicular (a type). Since the d l titanates were found to have a mirror plane 
perpendicular to the c axis, they have d type orbital order if they are JT distorted. 
In YTi0 3 , a small (3%) JT distortion was observed, but hardly any in LaTi0 3 . In 
their series of model Hartree-Fock calculations, Mizokawa and Fujimori could obtain 
an antiferromagnetic ground state for LaTi0 3 -but only by virtue of the spin-orbit 
coupling- as well as the ferromagnetic ground state for YTi0 3 . However, they could not 
explain why d- rather than a- type orbital order is realized in YTi0 3 [25]. This problem 
was later solved by Mizokawa, Khomskii, and Sawatzky by including A d orbitals 
together with the GdFe0 3 -type distortion in their model. However, these authors also 
noted that, given the smallness of any JT distortion in LaTi0 3 , this mechanism would 
make LaTi0 3 ferro- rather than antiferromagnetic [26]. 

Early on in the game, the correct magnetic orders had actually been computed 
by Solovyev, Hamada, and Terakura for LaTi0 3 [31] and by Sawada and Terakura for 
YTi0 3 [32] with the LDA+ U method. They adjusted U to the experimentally measured 
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Mott gap. Their orbital order for YTi0 3 was subsequently confirmed by NMR [33], 
neutron scattering [34], and resonant x-ray scattering [61]. 

The spin-wave spectra were then measured and found to be essentially isotropic and 
gapless in both titanates [10, 40]. This lead Khaliullin and collaborators to question the 
existence of orbital order in the titanates. They demonstrated that with the accepted 
orbital order in YTi0 3 , conventional super-exchange theory, as used in the model studies, 
leads to an anisotropic spin-wave spectrum and, even worse, to an antiferromagnetic 
ground state. If, on the other hand, the t 2g levels were degenerate, the isotropy of the 
spin-wave spectra could easily be explained by strong quantum fluctuations in orbital 
space (orbital liquid theory [16, 57]). The assumption that the t2 g levels are nearly 
degenerate was justified by the smallness of the observed JT distortions. 

Mochizuki and Imada [28] then pointed out that in order to obtain the observed G- 
type antiferromagnetism, as well as exchange-coupling constants in accord with the 
measured spin- wave spectra (J^ 1 = Jse°^ W ° = 15meV), it would suffice if the Ti 
3d 3z 2 _ x level were lower than the two other t 2g levels by an amount exceeding the 
spin-orbit splitting and kT, i.e. by >40 meV. A tiny deformation of the octahedron 
along [111] would enable this. They also gave arguments why such an orbital order had 
not been observed in resonant x-ray scattering [10]. At the same time [27], Mochizuki 
and Imada published strong-coupling (in U/t, as in super-exchange theory) studies of 
the antiferromagnetic-to-ferromagnetic transition, which takes place as the angle of the 
O-Ti-0 bond (~ degree of GdFeO"3-type distortion) decreases beyond that in LaTiO"3. In 
these studies, the JT distortion was assumed to be d-like, that is, with one level (xz on 
site 1) above a doubly-degenerate level (yz and xy on site 1). The dominant effect was 
found to be the t 2g -to-e g hopping via oxygen which increases with distortion, as we have 
seen in figures 14 and 15. Such hopping processes make it possible to reach virtual high- 
spin states with configuration t\ g e g and, hence, favour ferromagnetism. This study [27] 
therefore predicted that the orbital order is not influenced by the increasing GdFeC>3- 
type distortion, but that the exchange coupling along the z direction, where the t 2g -p-e g 
hopping is strongest, changes sign from A-type antiferromagnetic (i.e., ferromagnetic 
parallel to the mirror plane) to ferromagnetic. None of these studies [27, 28] included 
the A d orbitals which had been found crucial by Mizokawa, Khomskii, and Sawatzky 
[26]. 

Now, orbital-liquid theory had predicted that there be a significant contribution to 
the specific heat at low temperature, but this was not observed [18, 19]. Moreover, a 
3% JT stretch of one of the basal O squares into a rectangle was recently discovered 
in LaTiOs [12], thus bringing strong evidence for non-degenerate t 2g levels and a JT 
distortion of a different type than the one in YTi03. Mochizuki and Imada finally 
[24, 29] included the A d orbitals and realized that what stabilizes a single level 
with rf 3z 2 ii _ 1 -character, and therefore G-type antiferromagnetism, is not JT distortion, 
but the GdFe0 3 -type distortion which makes the A ions produce a crystal field. Its 
electrostatic and covalent contributions were estimated to give about the same lowest 
state, namely |1) = 0.68 \yz) + 0.41 \zx) + 0.61 \xy) . With this orbital order, they 
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obtained super-exchange coupling constants (J^? 1 = 19.7 meV and Jse°^ 10 ° = 18.5 meV) 
in almost perfect agreement with the experimental values [10]. Proceeding along the 
series of rare-earth titanates with decreasing ionic radius from LaTi03 to SmTi03, 
the J se values were found to decrease by 20%, and the lowest orbital to become 
|1) = 0.73 \yz) + 0.24 | zx) + 0.63 \xy) . From this continuous transition towards the 
case of YTi03, it was concluded that beyond Sm the JT distortion should be the factor 
which controls the magnetism, i.e. the JT distortion should be responsible for the 
ferromagnetism of YTi0 3 [29] . This work did not mention including the coupling to the 
e g degrees of freedom. 

Essentially the same conclusions were reached independently, and in one shot, by 
using a new, parameter-free density-functional approach [21]. The present paper tells 
the story as it is seen from that side. The fact that very similar conclusions concerning 
the splitting of the orbital degeneracy and its universal origin were reached in such a 
short time, not only by different theory groups, but by many experimental groups as well, 
strengthens the case for those conclusions. Nevertheless, we, too, must demonstrate that 
the orbital order obtained is consistent with the magnetic structures and the spin-wave 
spectra [10, 40]. We thus move a step backwards and repeat Khaliullin's calculation 
[40] of the inter- atomic exchange integrals by 2nd-order perturbation theory in t/U, 
this time using the crystal-fields and hopping integrals obtained from our LDA Wannier 
functions: 



6.2. Super- exchange couplings 

The super-exchange Hamiltonian is obtained from the t2 9 Hubbard Hamiltonian, 7~tt 2g , 
given by expression (27), including now the spin-flip terms neglected for the purpose 
of the DMFT calculations. Following the standard approach, we work in the many- 
electron representation which diagonalizes the on-site terms, Ti.^. R=RI + Ur, of Ti. t , 2g . 
The remaining, non-diagonal part of T~tt 2g is the LDA hopping, T. The super-exchange 
Hamiltonian, Tt se , is then obtained by Lowdin downfolding (A. 5) of the subspace which 
has one or more sites doubly occupied: 

n ae = V S T{(1 - V s ) {E - H t2g ) (1 - V^}' 1 TV S . (34) 

Here, V s projects onto the subspace of singly-occupied sites. To order t/U, the projection 
1 — V s can be substituted by Vd, which projects onto the space of one doubly-occupied 
site, and TC t2g can be substituted by the on-site part of the Hamiltonian. Expression 
(34) then reduces to the 2nd-order perturbation expansion for TC se used by Khaliullin 
et al. 

We now specialize to the cases of LaTi0 3 and YTi0 3 . As explained in the previous 
section, the main result of our many-body LDA+DMFT calculations was that the 
density matrix is diagonalized by the crystal- field orbitals, \m) , and that the orbital 
with the lowest energy, |1) , is occupied with 0.91 electron in LaTi0 3 and with 0.96 in 
YTi0 3 . In the following, we shall assume that |1) is completely occupied. For the two 
sites, R and R', between which we want to compute the super-exchange coupling, the 
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states with no double occupancy, i.e. t\ g + t\ , are then the Slater determinants |1) 
through 1 4) in the pictogram below, and the states with one double occupancy with 
which they couple, line by line in the pictogram, are the 15 t\ g states |5) through 1 19) 
(soft bras and kets denote two-electron states): 



|1) 



T I 



|5) 



|2) 



I T 
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|3) = 



T T 
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|18) = 
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|19) = 
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(35) 



In this space of 19 two-electron states, the matrix element of the Hamiltonian Ti. t2g 
between states (1| through (4| and |5) through 1 13) is a hopping integral, such as 
(lR\H t2g \7R') = (lR Hf° A 3R'^ = t 13 between crystal-field orbital |1) on site R 
and crystal-field orbital |3) on site Rl . The matrix elements in the £| g -subspace are the 
eigenfunctions of U plus the crystal- field splittings, A m = e m — e±. Now, the orders of 
magnitude produced by our LDA calculations for the titanates are: t ~ A ~ 0.3J ~ 
O.OAU, and we may therefore neglect A when evaluating the t\ g eigenfunctions (although 
this is not done in our numerical calculations). The t\ g eigenfunctions are then simply 
the p 2 -states listed following expression (27). The transformation from the two-electron 
basis (35) to the one which diagonalizes U is very simple: it merely involves combining 
states |6) and |8) into a singlet and a triplet (belonging to respectively 1 D and 3 P), 
the same for |7) and |9), and combining states |5), |14), and 1 15) into the 1 S singlet 
and the two singlets belonging to l D. One of the latter does not couple to the singly- 
occupied states |1) through |4), and the same holds for the three triplets formed from 
1 16) through |19). The singlet formed from states 1 16) and 1 17) belonging to 1 D does 
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not couple either. As a result, the off-diagonal block of H t2g — 2ei is then: 
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(36) 

Running along the rows, we have the following t^-states: the ^-state, three of the five 
^-states, and six of the nine 3 P-states. With the high-energy subspace diagonalized, 
it is now trivial to downfold to the low-energy subspace, |1) through |4), obtaining the 
super-exchange Hamiltonian (34) with E = 2e\. For the energy difference between the 
Neel state and the ferromagnetic state, we thus obtain: 



{3\Hse\3)-{l\H se \l) = t\ 



+ 



i 



12 



t 



13 



(37) 



1:1 \U + 2J ' i' - J J :l \ l J I - 3./ 

which for spin |, is | times the super-exchange coupling, J se , between the two sites 
considered. For simplicity, and because it is a good approximation, we have left out the 
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crystal-field splittings in the denominators, but they are trivial to include by returning 
to the full Hamiltonian (36). 

The super-exchange coupling (37) is the sum of terms caused by intra- and 
inter-level hops. As is well known, the first term is positive, that is, it favours 
antiferromagnetic alignment, and the latter is negative. Ferromagnetic alignment is 
thus favoured when the hopping to excited levels exceeds the hopping in the ground 
state by a critical ratio, r : 



For the actual value of the ratio between the Hund's-rule coupling to the Coulomb 
repulsion, J/U = 0.64/5 ~ 0.13, this equation gives r = 1.47. Similarly, for a given 
ratio of inter- to intra-level hopping integrals, one may find a critical value of the Hund's- 
rule coupling to be exceeded in order to produce ferromagnetic coupling: 



We now consider the super-exchange coupling between two nearest neighbours on 
either side of the mirror plane, i.e. those separated by 001, and between two nearest 
neighbours in the x|/-plane, i.e. those separated by 010 or 100. In the latter case, we 
must average over the two sites related by a glide mirror, and this amounts to taking the 
average of the 010 and 100 exchange couplings. Couplings between farther neighbours 
are very small because the squares of the corresponding hopping integrals determine 
J se , so we neglect them. In table 11, we have extracted t^m an d t\m f° r LaTiOs from 
equation (23), and t^m fr° m equation (25). For YTi03, we use equations (24) and (26) 
in the same way. From these hopping integrals and the values of U and J, we can 
now evaluate the various contributions in (37) to the super-exchange coupling, and the 
corresponding rows have been inserted in the tables. The summed up, total J se may 
be found in an inserted row. Finally, J se evaluated numerically without neglecting the 
crystal-field splittings has been included in a row marked J^ e . 

For LaTi03 and the coupling to 001, we see that the antiferromagnetic intra- 
level hopping dominates and that the total J®® 1 is 5.0 meV. Including the crystal- 
field splittings makes no difference in this case. Parallel to the mirror plane, a 
couple of large inter-level hopping integrals make the total exchange constant smaller, 
jioo/oio _ 2 2 me v, but still antiferromagnetic. These values are much smaller than the 
values obtained from the neutron scattering [10] and Raman experiments [70], and also 
than the values obtained with the model Hamiltonian [24, 29]. With the older structural 
data [44], we find the perpendicular and the parallel exchange-coupling constants to be 
both 6 meV. Using the low-temperature structural data [12], we compute the LDA 
Hamiltonian given in the left-hand side of table 7, and with that, we can calculate 
the crystal-field eigenvalues and vectors like in equation (19), and then transform the 
hopping integrals to the crystal-field representation like in equations (23) and (25). The 
result is given in table 12: due to a 10% increase of = the isotropy is increased, 




(38) 




(39) 
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Table 11. Hopping integrals, , and superexchange couplings, J se , in meV 



LaTi0 3 YTi0 3 



xyz 


001 


010 




100 




001 


010 




100 




1, m 


^lm Jse 


tlm 


Jse 


tlm Jse 




^lm Jse 


tlm 


Jse 


tlm Jse 




1,1 


105 9.1 


-98 


7.9 


-98 7.9 




-38 1.2 


-48 


1.9 


-48 1.9 




1,2 


31-0.2 


-192- 


7.0 


4-0.0 




-21 -0.1 


— 191 - 


7.0 


-84-1.3 




1,3 


143-3.9 


12 


0.0 


120-2.7 




97-1.8 


-130- 


3.2 


94-1.7 




Jse 


5.0 




0.9 


5.2 


3.0 


-0.7 




8.3 


-1.1 


-4.7 


° se 


5.0 








3.2 


-0.5 








-4.0 




0.99 








1.16 


1.85 








2.74 


(k 


0.19 








0.16 


0.10 








0.05 



Table 12. and J se in meV for LaTiOg at T=8 K. 
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1, m 
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J se 
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1,1 
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10.7 


-108 


9.6 


-108 


9.6 




1,2 


43 


-0.4 
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-6.3 


23 - 
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1,3 
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-3.8 
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-0.1 
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(k 
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Jse 1 = 6.5 and Jie ^ 010 = 5.0 meV, but our ab initio exchange coupling remains 2-3 
times smaller than the experimental value. Had we, instead of the computed crystal- 
field eigenvectors (19), used the simplified ones: 
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1 5 and 7 100/01 ° 



(40) 



10 meV, as 



where |1) = 3zf lx — 1, we would have obtained J\ 
shown in table 13 below. This demonstrates the extreme sensitivity of the exchange 
coupling constants to the crystal-field eigenvectors. 

For YTi03, the small intra-level 001 hopping makes the ferromagnetic inter-level 
coupling dominate so that J° e 01 = —0.5 meV, as seen in table 11. Parallel to the 
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Table 13. t\Z and J se in meV using the simplified eigenvectors (40) and (41). 
LaTi0 3 YTi0 3 



xyz 
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1, m 


t\m Jse 


t\m Jse 


tlm Jse 


tlm Jse 


tlm 


J se 


tlm 
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1,1 


150 18.5 


-128 13.5 


-128 13.5 


-35 1.0 


-84 


5.8 


-84 


5.8 




1,2 


51-0.5 


-133-3.4 


71-1.0 


-36 -0.2 


-132- 


3.3 


-100- 


1.9 




1,3 


126-3.0 


12-0.0 


101-1.9 


103-2.0 


-175- 


5.8 


65- 


0.8 




Jse 


15.0 


10.1 


10.610.3 


-1.3 




3.3 




3.1 


-0.1 




0.64 




0.71 


2.20 










1.49 


(k 


0.25 




0.24 


0.07 










0.13 



mirror-plane the large inter-orbital hopping integrals, however, make the calculated 
exchange coupling robustly ferromagnetic: j}e°^° W = — 4.0meV. This is consistent with 
the measured Curie temperature and the spin- wave spectra, which give —3 meV [40]. 
Had we, like for LaTi0 3 , used simplified crystal-field eigenvectors for YTi0 3 , namely 



(|1) |2> | 3 )) = (|y*> \xz) \xy)) 



/ 71 -71 ^ 
1 

v* 75 °y 



(41) 



we would have obtained J®® 1 = —1.3 and jjf '' 010 = —0.1 meV, as demonstrated in the 
right-hand side of table 13. So YTi0 3 does stay ferromagnetic, but just barely. 

From the hopping integrals, we can also calculate the ratio, r, between the intra- 
and inter-level hoppings using (38). The corresponding lines have been added at the 
bottom of tables 11-16. We see that for both titanates, r is nearly isotropic and, for 
LaTi0 3 , it is significantly smaller than the critical ratio, r =1.47, whereas for YTi0 3 , 
it is significantly larger. This agrees with the experimental facts. Obviously, the small 
tu and the large ti2 and ti 3 in YTi0 3 , make this material profit from the Hund's-rule 
coupling in the t\ g state and, hence, have a ferromagnetic ground state. 

Using (39) we can calculate the critical J/U ratio needed for the coupling to be 
ferromagnetic. (J/U) is seen to be significantly larger than the actual value (0.13) 
for LaTi0 3 , and significantly smaller in YTi0 3 . In the bottom panels of figure 24 we 
show the perpendicular (red) and parallel (blue) exchange couplings as functions of J/U, 
computed with inclusion of the crystal-field splittings. The latter are seen to have only 
minor effects, at least on the scale of the figure, but not on the scale of experiments, e.g. 
the largest deviation of the critical J/U values obtained without crystal- field splitting 
is found for the in-plane exchange coupling in LaTi0 3 . Here, the full calculation in the 
figure gives ( J/U) — 0.175, while the approximate value quoted in the table is 0.19. 
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Table 14. t\Z and J se in meV for cubic hopping. 



LaTi0 3 YTi0 3 
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tlm Jse 


tlm Jse 
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-20 -0.1 
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Table 15. t\Z and J se in meV for simplified eigenvectors and cubic hopping. 



LaTi0 3 YTi0 3 
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1, m 
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Jse 
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1,1 


100 8.2 


-100 8.2 
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75 4.6 


-75 


4.6 


-75 4.6 


1,2 


0.0 


-61-0.7 
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1.1 


-75-1.1 


1,3 


71-1.0 


35-0.2 


35-0.2 


0.0 
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2.1 


0.0 


J. se 


7.3 


7.3 


7.3 7.3 


3.6 




1.4 


3.6 2.5 




0.50 




0.50 
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1.00 


(k 


0.28 




0.28 


0.24 
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Our LDA-NMTO results should be seen on the background of the claim [40] that 
orbital order is inconsistent with the observed magnetic ground states and isotropic spin- 
wave spectra in YTi03. This claim was based on the assumption that the hopping is as 
in the cubic structure, where the xy, yz, and xz bands are identical, independent, and 
2D. In that case, there is only one nearest-neighbour hopping integral, the effective t^ 
between two similar orbitals in the same plane (in case of SrVC>3, this is the —281 meV 
quoted in table 1). As we have seen and explained in great detail in section 3, the 
hopping between the xy, yz, and xz LDA Wannier functions is, however, very different 
due to the GdFe0 3 -type distortion in the titanates. In order to demonstrate that this 
is the reason for the discrepancy between our calculation and that of Khaliullin [40], 
we first of all note that our expression (37) reduces to their equation (3) if one assumes 
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Figure 24. Nearest-neighbour exchange coupling constants perpendicular (red) and 
parallel (blue) to the mirror planes as functions of the ratio between the Hund's-rule 
coupling and the on-site repulsion, J/U. See figure 4 in reference [40]. Antiferro- and 
ferromagnetic couplings are respectively positive and negative. Input were U = 5 eV 
and the orbitally ordered states given by equations (19) and (20). In the upper row, 
degenerate ^g-levels and cubic hopping with t- K =— 150 meV were used, like in table 
14. The values of the exchange-coupling constants for J/U — 0.64/5 = 0.13 are those 
listed in table 11 as J^,. 

cubic hopping and uses the eigenvector of the lowest crystal-field state |1) in YTiOs. 
Since equation (37) is general, we can use it also for "cubic" LaTiOs. In table 14 
we have therefore listed the hopping integrals obtained by applying our crystal-field 
transformations (19) and (20) to cubic hopping with t n = —150 meV. In this case, the 
difference between LaTi03 and YTiOa stems exclusively from the different crystal-field 
eigenvectors. The differences with the ab initio values of the hopping integrals listed 
in table 11 are seen to be large, although some of the trends are captured by the ab 
initio crystal-field eigenvectors. From the values of r and (J/U) , we realize that, 
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although both "cubic materials" are antiferromagnetic, as claimed in [40], YTi0 3 does 
have a tendency to couple ferromagnetically in the plane. The calculation for the cubic 
titanates is included as the top panel of figure 24. 

If we finally take cubic hoppings, as well as the simplified eigenvector (40) for 
LaTi03 and (41) for YTiOs, these then being the only difference between the two 
materials, we obtain the results shown in table 15. Both "materials" are fairly isotropic, 
G-type antiferromagnets with the exchange coupling two times larger in LaTi0 3 than 
in YTi0 3 . 

In conclusion, the sign of the super-exchange coupling that we calculate using 
conventional theory and LDA t2 9 Wannier functions is consistent with the observed 
magnetic ground states of LaTiOs and YTiC>3. The mechanism by which the different 
magnetic orders come out of the calculations is through orbital order caused by the 
Coulomb repulsion and a strong crystal-field, and through distortion of the cubic 
hopping integrals. Both the crystal field and the modification of the hopping integrals 
are caused by the GdFe0 3 -type distortion. The fact that we can explain both the 
Mott gap and the observed magnetic order from the orbital order and hopping integrals 
calculated without adjustable parameters with the LDA, seems to rule out the orbital- 
liquid scenario for the t\ g titanates. 

The discrepancies between experiment and computation do deserve further 
investigation, however. 

The size of our J^? 1 exchange-coupling for YTiOs agrees with the experimental 

values, but the anisotropy is too large. Due to the near-cancellation of terms, the 

calculated exchange coupling constants are very small, as they should be in order 

to account for the low Curie temperature, but this also makes them sensitive to 

detail. For instance, the values of our calculated 2nd-nearest-neighbour constants, 
jon/ioi = +Q 2 meV; joii/ioi = +Q 3 meV) juo = +Q 5 me y 5 and jiio = +Q 2 meV; are 

antiferromagnetic and not negligible compared with the lst-nearest-neighbour values in 
table 11. We shall see a further example of this sensitivity in section 6.3, where we 
consider the influence of the JT distortion. Moreover, the computational inaccuracies 
mentioned in appendix B may influence these delicate exchange-coupling constants. 

For LaTiOs, our exchange couplings in table 12 are isotropic, but three times 
smaller than the experimental values. In this case, our calculated 2nd-nearest-neighbour 
constants, J s ° e 11/101 = -0.2 meV, J° e 1I/101 = -0.1 meV, J™ = +0.9 meV, and 
•^ie 10 = meV, are small in comparison with the lst-nearest-neighbour values. Also 
the computational inaccuracies mentioned in appendix B are presumably too small to 
explain this discrepancy. 

In general, one suspects that the LDA overestimates covalency. If that would cause 
the scale of our hopping integrals to be too large, our value of U chosen to account for 
the Mott transition would also be too large in order to compensate for this, but then 
our t 2 /U should be too large - but, in fact, the opposite seems to be true for LaTi0 3 . 
Nevertheless, since our hopping integrals result from a delicate balance between Ti-O, 
Ti-Ti, A-O, and Ti-A covalencies as seen in sections 3.3, 3.4, and 3.7, the degree to which 
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the LDA gets this balance right is decisive for the accuracy of the calculated exchange 
couplings. To check this further, calculations with other density functionals are needed, 
and so are experiments to test the implications of the strong cation covalency predicted 
by our [21] and earlier model calculations [24, 26]. 

Another reason for the factor-three discrepancy in LaTiO,3 could be that the e g 
degrees of freedom should have been treated explicitly because, as was pointed out by 
Mochizuki and Imada [27], the energy of the e g band is merely a few times 2 J above 
that of the t 2g band (see figure 4). However, not fully treating the e g degrees of freedom 
should lead us to underestimate the tendency towards ferromagnetism - a problem that 
we do not have for LaTiOs. Besides, if we were to include the e g degrees of freedom 
in the Hubbard model, we could do so by adding the space of the Wannier functions 
for the e g band to the space of Wannier functions for the ti g band, because with this 
choice the Wannier orbitals and, hence, the parameter values for the t2 g band would be 
unchanged. But in this d basis there is no single-particle coupling between the t 2g and 
e g subspaces, and therefore no coupling in the 2nd-order perturbation expression (34). 
As a result, our super-exchange calculation in which only the t 2g degrees of freedom 
are treated and the e g - as well as all other single-particle channels are downfolded, is 
correct. 

The factor-three discrepancy might be related to the fact that the magnetic moment 
in LaTiOs is anomalously small (0.57 hb)- Our assumption of complete orbital order and 
subsequent use of 2nd-order perturbation theory are clearly insufficient if the orbitally- 
ordered state changes at low temperature. A future spin-polarized LDA+DMFT 
calculation should allow us to calculate the magnetic moment and to keep full account 
of the inter-site coupling between spin and orbital degrees of freedom. For a general 
discussion of the why the magnetic moment is so small in LaTi0 3 , we refer to reference 
[29]. 

Finally, it is possible that in these ti g materials with strong A covalency the 
eigenvectors of the density matrix calculated more accurately than in the single-site 
DMFT differ sufficiently from those of the crystal-field Hamiltonian to influence the 
exchange-coupling constants. 

6.3. Influence of the JT distortion in YTiO^ 

We have found that the JT distortion in YTiC>3 does not create, but merely reflects 
the orbital order driven by the Coulomb correlations and controlled by the GdFe03- 
type distortion through the crystal field set up by the A ions. Mochizuki and Imada 
[29] now speculated that for strongly GdFe0 3 -distorted materials like YTi0 3 , where the 
magnetism is weak, the JT distortion could be a controlling factor for the magnetic 
order. This we now investigate. 

Using the LDA Hamiltonian (table 7) computed for our hypothetical YTi0 3 without 
JT distortion, we obtain the hopping integrals and exchange couplings given in table 16: 
Jfg 1 = —3.0 meV and Jle°^ 0W = +2.5 meV. This gives ferromagnetic order along z and 
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antiferromagnetic order in the plane, that is, C-type antiferromagnetism and, hence, a 
completely different result than obtained for the proper structure. 

For the high-pressure phase observed by Loa et al. to have strongly reduced JT 
distortion, but slightly increased GdFeC-3-type distortion, we find the results given in 
the right-hand side of table 16: J° e 01 = +4.9 meV and J S T /01 ° = -7.8 meV. This 
causes a third kind of ground state with robust antiferromagnetic order along z and 
ferromagnetic order in the plane: a type- A antiferromagnet. 

In conclusion, our ab initio calculations support the speculation of Mochizuki and 
Imada [29] that, although not very important for the orbital order, the JT distortion 
is of crucial importance for the magnetic order in YTiC>3. This should be tested 
experimentally. 

7. Unfolding the orthorhombic band 

In this paper we have first seen three degenerate, independent, and nearly two- 
dimensional cubic bands develop into twelve inequivalent and coupled orthorhombic 
bands, which narrow down and develop a pseudo-gap under increasing GdFeOs-type 
distortion (figures 3, 10, 17, and 18). In order to treat the Coulomb correlations, 
we have constructed a set of highly localized t2 g Wannier functions for these bands 
and have defined the corresponding Hubbard Hamiltonian. With this set of Wannier 
functions, the splitting of the lower t 2g levels reaches a mere 10% of the bandwidth, 
but this, together with the 30% decrease of the entire bandwidth through the series 
(table 8) is important for determining where in the series the Mott transition occurs 
(table 10). However, U' c /W decreases by 50% when going from SrVO"3 to YTiO"3, and 
that is significantly more that the ~30% which can be gained by a crystal-field induced 
decrease of the effective degeneracy from 3 to 1. Therefore, the width of the lowest 
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subband must decrease faster than the width of the entire t 2g band. 

A second surprising result of the LDA+DMFT calculations was that once U' exceeds 
~ W, the eigenvectors of the density matrix are essentially the same as those of the 
on-site LDA Hamiltonian, and as U' ~ U' c , the orbital fluctuations become strongly 
suppressed. As a result, only the orbital of lowest energy is occupied in LaTi0 3 and 
YTi0 3 . As seen in figure 21, the (1 |e| 1) element of the spectral function is divided 
into a lower occupied and a higher unoccupied Hubbard band, while the diagonal 
elements of the higher orbitals, (2 \e\ 2) and (3 |e| 3) , have weight almost exclusively 
in the upper, unoccupied Hubbard band. The off-diagonal elements, (1 \e\ 2) , (2 |e| 3) , 
and (1 |e| 3) , are completely negligible. We thus found orbital order, particularly in 
YTi0 3 . Our conclusion of those two findings was that the self-energy must behave in 
the way explained around equations (31) and (32). 

In section 6 we found that the increasing tendency towards ferromagnetic coupling 
is due to the increasing strength of the hopping in the orbitally ordered state from the 
lowest crystal-field orbital to the higher orbitals on the neighbouring sites, compared 
with that of the hopping between the lowest orbitals. This is what the ratio r given in 
tables 11-16 measures. Increasing r, of course, tends to decrease the width of band I. 

What remains unexplained is the 50% effect brought out by table 10 or, in other 
words, why for merely U' > 1.6W in LaTi0 3 and 1.3VF in YTi0 3 the lower half of 
one band (per Ti) - considerably more narrow than W - lies below the bottom of 
the others, and why this band has pure LDA crystal- field character |1) . In figure 17 
it is hard to identify such four orthorhombic LDA bands with |1) character. That 
is, the orthorhombic LDA bands display little tendency towards symmetry breaking. 
Nevertheless, we shall now demonstrate that behind the LDA pseudo-gap in the 
orthorhombic GdFe0 3 -distorted perovskites there is a real gap. This gap is a direct 
one in a pseudo-cubic k space where it splits off the lowest t 2g band. 

The way we arrived at this result was, first, to see whether with the NMTO method 
we could produce one Wannier function which describes the four lowest orthorhombic 
bands and has the character of the lowest crystal-field orbital. The result of this attempt 
is demonstrated in figure 25. The dashed black bands are the correct LDA bands, 
as obtained with a large NMTO basis set, and the 4 red bands are obtained with a 
truly minimal NMTO basis set containing only one t 2g orbital per site, the one with 
the symmetry of the lowest crystal-field orbital, i.e. with the eigenvector |1) defined in 
equations (19) and (20) for LaTi0 3 and YTi0 3 , respectively. Since this orbital is not the 
partner of two other t 2g orbitals, which together with it would span the entire t 2g band, 
the 1 2) and |3) characters (figure 14) are downfolded into the tail of this new, extended 
|1) orbital, which we shall name |I) and show in figure 26. Moreover, the energy mesh of 
1 1) was chosen to span only the occupied part of the t 2g band. The blue bands in figure 
25 are obtained from another truly minimal NMTO basis set with two t 2g orbitals per 
site, | II) and |III), with respectively |2) and |3) character, and with |1) downfolded in 
the tails. Its energy mesh is mostly in the empty part of the t 2g band (appendix A). 
The marvel is that orbital |I) picks 4 of the 12 t 2g bands, that orbitals |II) and | III) 
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Table 17. H LDA in meV for LaTi0 3 [12]. 
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Table 18. H LDA in meV for YTi0 3 . 
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Figure 25. LDA bands obtained with a large NMTO basis set (dashed black) and 
with truly minimal sets, downfolded to respectively the lowest ti g crystal-field orbital 
1 1) (red) and the two highest crystal-field orbitals |II) and | III) (blue). The following 
NMTO energy meshes were used: for the LaTiC>3 red bands e„ =—0.55, —0.45, and 0.04 
eV, and for the blue bands e n =-0.94, 0.70, 1.11 and 1.17 eV. For YTi0 3 , e n =-0.25, 
-0.12, and 0.40 eV for the red bands, and e„ =-0.80, 0.70, 1.04, 1.31, and 1.65 eV for 
the blue bands. 



together pick the remaining 8 t 2g bands, and that the 4th band near T is the lowest of 
the blue bands and has exclusively |III) character. The fact that such a decomposition 
of the bandstructure is possible is not trivial, at least not from the point of view of the 
orthorhombic bands. Since all the LDA Bloch states are orthonormal, so are the three 
extended crystal-field Wannier functions, to the extent that they span the LDA bands, 
and -as figure 25 demonstrates- they do this with good accuracy. 

Comparison of the bands in figure 25 with the ones decorated with the local 1,2,3 
characters in figure 17 confirms that band I picks four low bands with dominant 1- 
character. The orbital which describes those four bands has 1-character in its head and 
the minority characters 2 and 3 in its tail, on the neighbouring Ti sites. The way in 
which it acquires a low energy is for a given neighbour to pick that linear combination 
of orbitals 2 and 3 which provides the same orientation of the oxygen p function to 
be shared, and then to add this linear combination with the opposite phase and an 
appropriate amplitude. In this way, the resulting oxygen p function will bond with the 
Ti characters on the neighbours. 

The Hamiltonian for band I and the one for bands II and III are now 
orthonormalized and Fourier transformed. The resulting on-site and hopping matrix 
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elements are given in tables 17 and 18. Remember that, by virtue of our NMTO 
construction, there is no LDA interaction between orbital I and the two other orbitals. 
We note that the splitting between level I and the nearly degenerate levels II and III 
is huge, 0.7 eV. The bands have therefore really been separated and what we have is 
extreme orbital order at the LDA level. For LaTi0 3 (YTiC^), the energy of orbital I 
has been lowered by 0.35 (0.30) eV from that of orbital 1, and the energies of orbitals 
II and III have been raised from those of orbitals 2 and 3 by, in total, the same 
amount. This has been achieved by adding to orbital 1 bonding 2- and 3-characters 
at the neighbouring sites, as was mentioned above, and by adding to orbitals 2 and 3 
antibonding 1-character at the neighbours. It is therefore also obvious that we have got 
rid of inter-orbital hopping at the expense of inter-site Coulomb repulsion and Hund's- 
rule coupling, without which the ferromagnetism of YTi0 3 cannot be explained. For the 
set of extended crystal-field orbitals, the Fermi level almost coincides with the lowest 
level; it is slightly below in LaTi0 3 and a bit above in YTi0 3 . This is because the 4th 
band at T (band III) dips deeper below the Fermi level in LaTi0 3 than in YTi0 3 . 

One might argue that these Wannier orbitals are far too extended to be of any use 
for describing correlated electrons. The way to think about these orbitals is, however, 
as the initial, S=0 orbitals in a self-consistent LDA+DMFT approach like the one 
recently suggested [51]. As self-consistency is approached, the self-energy will separate 
the bands and the extended Wannier orbitals will localize so that |I) — > |1) , |II) — > |2) , 
and | III) — > 1 3). So what the initial orbitals do is to tell us in which way, if any, "the 
chemistry" wants to break the symmetry. In the present case of t\ g titanates, it is 
clear that the chemistry tells us that we have one half-filled band, rather than three 
degenerate |-full bands. For this one band the DMFT will finally create a self-energy 
which will separate it into a lower and an upper Hubbard band, provided that the width 
of this band is sufficiently small. 

For orbital I, the lst-nearest neighbour hopping integral perpendicular to the mirror 
plane is positive, while those parallel to the plane are negative, and merely half the size. 
Moreover, the hopping integrals in YTi0 3 are half the size of those in LaTi0 3 . A way 
to analyze the hoppings would be to express |I) as a linear combination of |1) , |2) , and 
1 3) on the same site and on the 6 neighbour sites. This would, for instance, explain why 
the hopping is twice as strong in the z direction as in the parallel direction, despite the 
fact -obvious from figure 26- that the orbital is most extended in the y direction and 
least in the z direction. As we have experienced before, hopping integrals are sums of 
terms with alternating signs. 

For LaTi0 3 and YTi0 3 we have thus succeeded in decoupling one orbital and in 
removing the coupling between the two other orbitals approximately. The consequence 
of this is that we can fold the band structure out in a pseudo-cubic Brillouin zone (see 
figures 2 and 9). The reason is that the group of covering operations for a single orbital 
on the B sites is cyclic, provided that the orbital is defined as in the bottom row of figure 
8. As a result, we can now let expik ■ R be the irreducible representations of this cyclic 
group with R being primitive monoclinic translations. The meaning of k in the large, 
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Figure 26. Extended crystal-field NMTO |I) spanning the red LDA bands of figure 
25 and with the same local character as |1) in figure 14Compared with that figure, the 
present viewpoint is more from the top, looking down on a mirror-plane. The central 
part of the orbital is in subcell 1 000 and the rest is in subcells 2: 010 to the left 
and right, and 100 up and down. The amplitude in subcell 3 001 is smaller and has 
been truncated for the sake of clarity. The contour is the same as in figure 14. This 
orbital (before orthogonalization) has 63% of its charge density on the central AgBOg 
unit, and this holds for both titanates. When moved to the neighbouring Ti sites, this 
orbital follows the space-group symmetry like orbital |1) in figure 16 (orbital order). 



primitive monoclinic BZ is thus different from that in section 3.6, where we combined 
the B-centred orbitals, |R, m) , into Bloch sums following the middle-row convention of 
figure 8. In the present, so-called pseudo-cubic scheme, we still use equations (8) and 
(9), but there is only a single orbital shape, |R) , and there will be no coupling between 
the four Q vectors. 

For band I of LaTi03 and YTiC>3 the unfolding is demonstrated in figure 27. For a 
single orbital, the band dispersion is given by equation (1), an amazingly simple result 
considering the complicated orthorhombic bandstructures at the bottom of figure 27, 
not to speak of figure 10. Since the dominating nearest-neighbour hopping integral, £ 001 , 
is positive and the perpendicular ones, t 010 =t 100 , are negative, the bottom of the band 
is at 007r and the top is in the k z =0 plane near 7T7tO and 7r00. This may seem unfamiliar, 
but is related to the necessity of using the orbital-order convention at the bottom row 
of figure 8 instead of the physical convention in the middle row. We shall return to this 
later when we consider the cubic t 2g bands 

In figure 28 we include the hybridized bands II and III, this time in the orthorhombic 
zone because the hybridization cannot be folded out. The result is the dashed black 
band structure, which is the t% g band structure, but with all hops longer than to the 
120 neighbours truncated. The red band structure is that of band I, which is identical 
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to the four lowest black bands, and to the bands in the bottom panel of figure 27. The 
green and blue bands are, respectively, the unhybridized bands II and III, and they are 
seen to be reasonably accurate below, and up to 0.3 eV above, the Fermi level. 

With the hybridization between bands II and III neglected, also these bands may 
be folded out: all three bands are then given by equation (1), and they are shown in 
figure 29. We see that the lowest band is separated from the two others by a direct gap. 
Hence, the pseudo-gap is a real gap in pseudo-cubic k space, but with the lowest and 
the upper bands overlapping in energy. 

Before continuing the discussion, let us first explain how the pseudo-gap arises under 
increasing GdFeC>3-type distortion. This we did in section 3.6 by using conventional Q 
coupling and considering the orthorhombic S G point. Now, we want to use the pseudo- 
cubic k representation and connect back to the cubic limit, with its three degenerate, 
non-interacting, two-dimensional bands. That limit is given in the bottom row of figure 
29, where we used the 2nd-nearest-neighbour model (2) with t$ = 0. However, the 
pseudo-cubic representation is valid only if we neglect hybridization between the three 
Bloch waves defined with orbital order (bottom row of figure 8). This means that we 
first have to break the cubic symmetry and prepare the cubic bands for the GdFeC>3- 
type distortion. In other words, we have to pick band I, and this will of course be the 
xy band because this is the band for which the orbital order does not break the physical 
coupling (for xy, the bottom and middle-row definitions in figure 8 are the same). The 
k x , k y dispersion is thus given by equation (2), also in the pseudo-cubic zone; this is the 
red band in figure 29. With t$ = 0, this band has no k z dispersion and its minimum is 
at 00k 2 , the saddle-points at n0k z and 0irk z , and the maximum at irnk z . The bandwidth 
is 8 l^l , and with increasing f a /t n = r, the energy of the saddle-points shifts away from 
the minimum towards the maximum, i.e. the lower part of the band is stretched. The 
value of r chosen for the figure is that of SrV0 3 , whose cubic band structure was shown 
along the same path in figure 11. Moreover, t n = —250 meV, which is numerically a bit 
smaller than that of SrVC>3, but is like those of CaVO^. 

Next, we consider the yz and xz bands. Lets us be pedagogical and perform their 
unfolding explicitly. The hopping Hamiltonian in the four-site orthorhombic Bloch 
representation is given by: 
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Figure 27. Folding band I out from the orthorhombic (bottom) to the pseudo-cubic 
BZ (top). See figure 9. In the orthorhombic zone with 4 sites, R, the eigenvectors are 
^expiQ • R with Q = 000 (red), 007T (blue), 7T7tO (green), and tttttt (black). The unit 
along the abscissa is n. The bands were produced using tables 17 and 18. The pseudo- 
cubic bands are given by equation (1). The letters denote high-symmetry points in the 
orthorhombic BZ and the path taken along the simple monoclinic BZ is the same as 
in figure 11, but the pseudo-cubic k includes the orbital order. 



where the rows and columns have been ordered so that the first 4 Bloch waves belong to 
the middle, and the last 4 Bloch waves to the right column in the bottom row of figure 
8. After the unitary transformation (9), H/ (2 |t„.|) becomes: 
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where for simplicity of notation we have defined 



u 



1 / \ , 1 

- {c x + Cy) , and r - {< - < tl 



2 v_ * ^ y> ' 

Along the diagonal we only have the k z dispersion, as is also obvious from the bottom 
row of figure 8. We may, however, get rid of the off-diagonal elements proportional to 
u by transformation to 

|III) = -L (\ yz ) - \ xz )) ail d |II) = -^(|yz) + !**», 

where \yz) and \xz) now refer to subcell 1, or, equivalently, follow the notation of the 
bottom row in figure 8. This yields a block-diagonal Hamiltonian with the following 
four blocks: 
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where according to equation (2), 



e(k) 



— 2 \- e v z 



k)+e a 



(k)\ = 2t n (u + c z ) + U' a uc z , 
(k)} = -2t n (l + 2rc z )v, 



are, respectively, the average and half the difference of the yz and xz band dispersions 
with k in the orthorhombic zone. The eight yz and xz bands are thus obtained by 
hybridization between pairs of degenerate average bands, II (Q, k) and III (Q + ttttO, k). 
There exists such a pair of hybridized bands, e yz (k — Q) and e xz (k — Q) , for each of 
the four Q vectors. These are 8 of the 12 bands seen in figure 11 for SrV0 3 , and the 
remaining 4 are the e xy (k + Q) bands. Note that the e XiXj (k) and e XiXj (k — Q) bands 
are degenerate at the corresponding boundary of the orthorhombic BZ (red in figure 9). 

The orthorhombic zone may now be folded out in the k z direction because the 1st 
and 2nd Hamiltonians (43) are identical with, respectively, the 1st and 2nd Hamiltonians 
(42), if considered as functions of the wave vector k = k + Q 2 . If we therefore let the 
latter run over the orthorhombic zone doubled in the k z direction, the Hamiltonians (43) 
should be dropped. This is exact for any band structure, and examination of any set 
of orthorhombic bands in this paper will reveal that they can be folded out along any 
vertical path, i. e. a pair of bands along T Z can be folded out to one band along 000 



007r, and similarly for Y T 



2 2 







7T 7T 

2 2 



it, and S R — > 7r00 - ttOtt. This is also the 



reason why all orthorhombic bands are at least 2-fold degenerate on the horizontal face 
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Figure 28. Hybridized (dashed black) and unhybridized I (red), II (green), and III 
(blue) bands produced from the on-site and hopping integrals in tables 17 and 18. 



k z = |. The cubic bands folded out to the double orthorhombic zone are the yz and xz 
bands (2) translated to tttttt and the same bands translated to 007r, as can be seen from 
the Hamiltonians (42). 

Folding out to a pseudo-cubic BZ is not exact for a Hamiltonian, but only for 
a single band. Consider, e.g., the cubic band Hamiltonians (42) as functions of the 
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We have obviously succeeded in folding out the diagonal elements, but the hybridization 
cannot be folded out. In the pseudo-cubic zone, we have two bands, namely the average 
band, e (k) , translated to the sites 007r and tttttt. But there is no way in which these two 
bands can hybridize in the pseudo-cubic zone to yield the proper cubic bands. That can 
only occur after translation by ttttO, i.e. by returning to the double orthorhombic zone. 
Although important for nearly cubic bands, the neglected coupling between bands II 
and III causes no problem for the lowest | of the t^g band in the titanates. This must be 
so because we have succeeded in downfolding band I correctly with the NMTO method. 

The cubic bands are given in the bottom part of figure 29. In red e xy (k), in green 
e yz (k — 7T7T7r) and e xz (k — tttttt) , and in blue e yz (k — 007r) and e xz (k — 007r). These yz 
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and xz bands are degenerate, and therefore equal their respective averages, e (k— 7T7T7r) 
with wave function yz + xz and e (k— 007r) with wave function yz — xz, in the (±110) 
planes containing 007r, 000, ±7T7rO, and ±.-mm. 

In this pseudo-cubic representation we may now follow the development of the bands 
as a function of the GdFe0 3 -type distortion. For this purpose it is helpful to use the 
orthorhombic labelling at the bottom of figure 29 to relate the simple pseudo-cubic bands 
to the cubically averaged bands in figure 11 and to the projected orthorhombic bands 
in figure 17. The pseudo-cubic bands of CaV0 3 are intermediate between those of the 
cubic 2nd-nearest neighbour model and those of LaTi0 3 . Starting from the cubic bands, 
the 45meV coupling between the yz and xz orbitals in CaV0 3 (table 2), produces a 90 
meV gap between the green and blue bands at X Q | ^tt, and similarly at Y Q f f vr, just 
above the Fermi level. The couplings between the other pairs of orbitals are ~30 meV 
and the 3-fold degenerate red-green-blue level at T G ||| splits well above Sp. Along 
T 007r - S G nOn - R G 7r0| - S G 7r00 there is a splitting between the red xy and the 
green band almost at the Fermi level, and at high energy there is a splitting of all three 
bands. Moreover, the :q/-band is lowered by ~70 meV due to the crystal-field splitting. 
Finally, going along r o -rtttt - Z D -kit | - T a ttttO - X ||0 - r o 000 there is a red-blue, a 
green-blue, and a red-green splitting. Of these, the latter splitting between the xy and 
the yz + xz bands is at the lowest energy. In conclusion, due to the coupling between 
the xy, yz, and xz orbitals in CaV0 3 , there are small splittings (<100 meV) at the 
crossings of the cubic bands, and this separates off a lowest band. The top of this band 
is a red-green maximum between X Q ||0 and T Q 000, and the corresponding width of 
the lowest band can be read from figure 10 or 17 as being 1.80 eV. This is substantially 
smaller than the t 2g bandwidth (2.45 eV) listed in table 8. Due to the sharp avoided 
crossings, the Wannier function for the lowest band in CaV0 3 has very long range and 
is certainly not a suitable basis for a Hubbard model. 

We now go to the titanates, for which the red, green, and blue bands in figure 29 are 
the I, II, and III bands corresponding to the extended Wannier functions defined earlier 
in this section. The strongest nearest-neighbour hopping integral coupling between the 
xy, yz and xz orbitals (tables 4 and 5) has increased to 75 meV in LaTi0 3 and 83 meV 
in YTi0 3 , and the avoided crossings at X G f fvr, Y Q |f7r, and T ||| have increased 
correspondingly to 150 and 166 meV. This can be seen in figures 10 or 17, but not really 
in figure 29 where the truncation of the separate Fourier series (1) for the red and green 
bands has rounded off the avoided crossings. Moreover, as may be seen from figure 
25, numerical truncations in the NMTO calculation of the extended Wannier functions 
also produce some error in the II+III Hamiltonian near sharp avoided crossings, most 
noticeably at X Q . Along T Q 000 - Z D 00^ - r o 007r there is now a red-blue avoided 
crossing which was not noticed in the 2nd-nearest-neighbour model where the xy band 
neither disperses along k z , nor couples to the yz and xz bands. As seen from figure 17, 
this avoided crossing is between bands I and III near T OOn. Here, band III is purely 3- 
like and much more dispersive in LaTi0 3 than in YTi0 3 , and it also lies higher because 
the 1-3 crystal field splitting is 205 meV in LaTi0 3 , but 330 meV in YTi0 3 . These 
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Figure 29. Bandstructures in the pseudo-cubic BZ (figure 9). For the titanates, the 
red band (I) is the same as in figure 27, and the green (II) and blue (III) bands are 
identical with those in figure 28, but folded out. The hybridization between II and III 
is neglected. The bands are described by equation (1) with the parameters given in 
tables 17 and 18. For the cubic model, given by equation (2) with t v = —250 meV, 
t'^/t^ = r = 0.34 and tg = 0, the red band (I) is xy, the green bands are e yz (k— tttttt^ 
and e xz (k— tttttt^ , and the blue bands are e yz (k— 007r) and e xz (k— 007r). These yz 
and xz bands are degenerate in the (±110) planes containing 007r, 000, ±7T7tO, and 
zLtttttt. Here, they therefore equal their respective averages, namely e (k— tttttt^J , the 
green yz + xz band II, and e (k— 007r) , the blue yz — xz band III. The cubic bands are 
roughly those of SrVC>3 and are shown along the same path in reciprocal space as in 
figure 11. The unit along the abscissa is 7r. 
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Table 19. Width of the lowest band, I, in eV and the ratio Wi/W 



Wi 

Wi/W 



SrV0 3 [42] CaV0 3 [43] LaTi0 3 [44] LaTi0 3 [12] YTiO 3 [20] 
2.85 1.80 1.40 1.30 0.95 

1.00 0.74 0.63 0.62 0.50 



differences in dispersion and position of band 3 might be the underlying cause for the 
differences in the onset of optical absorption calculated in subsection 5.5 and discussed 
most recently by Ruckamp et al. [63]. 

Avoided crossings give rise to small tongues in the joint density of states, 



extending down to ~0.15 eV, but the strong onset of Jf (u) occurs at higher energies: 
at 0.24 eV in LaTi0 3 , and at 0.33 eV in YTi0 3 . This relatively large gap is due to the 
combined effects of orbital order, strong hopping from orbital 1 to orbitals 2 and 3, and 
comparatively weak hopping between orbitals 1 (table 11). This gives band I a relatively 
small width, and it causes orbital I to have a much lower energy (0.7 eV) than orbitals 
II and II. This energy is the number of nearest neighbours times twice a typical hopping 
integral between orbital 1 and orbital 2 or 3, times some reduction factor due to the 
localization of the extended Wannier functions. This may be realized from figure 26 
and the discussion given at the beginning of this section. Even in the basis of xy, yz, xz 
orbitals, the inter-orbital hopping increases strongly along the series, a fact mentioned in 
connection with table 6. Going now the crystal-field basis dramatically enhances inter- 
orbital hopping. Specifically, the ratio r defined in equation (38) is larger for YTi0 3 
than for LaTi0 3 (table 11). This is what makes YTi0 3 ferromagnetic and what tends 
to make the relative width of band I small. 

The pseudo-gap penetrates the density of states in figure 18 better in YTi0 3 than 
in LaTi0 3 , because there is less indirect band overlap in YTi0 3 . The chemical reason 
is that the maxima of band I between T Q ttttO and X ||0, near S G 7r00, and at T G ||| 
are pushed down by the Ad 3z 2 _ x and d xy characters, as was demonstrated in figure 12. 
In terms of effective hopping integrals, this is what makes in YTi0 3 anomalously 
small, as was explained in connection with equation (5) and referred to many times. 
This then points to the common reason why YTi0 3 is ferromagnetic, both according 
to super-exchange theory and in the Stoner band picture. Finally we can read off the 
width, Wi, of the lowest band in figures 10 and 17, and list them in table 19 together 
with the relative subband widths W\/W . 
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8. Summary and Outlook 

By means of Wannier functions obtained from ah initio density-functional (LDA) 
calculations we have studied the series of orthorhombic perovskites AB0 3 = SrV0 3 , 
CaV0 3 , LaTi0 3 ,and YTi0 3 in which, under the influence of an increasing GdFe0 3 -type 
distortion, the single B t 2g electron becomes increasingly localized and undergoes a Mott 
transition between CaV0 3 and LaTi0 3 . The energy bands in figure 1 and the Wannier 
functions for the oxygen 2p band in figures 5 and 6 show that covalency between the 
occupied O p and the empty large-cation A s and, in particular, d states is an important 
mechanism of this distortion (figure 7). The A d states pda-bond with the same oxygen 
p orbitals as those with which the empty B t 2g states pc?7r-bond and, as a result, the B 
t 2g Wannier functions have, not only oxygen pdn but also residual A d 3z 2^_ l and d xy 
characters (figures 8, 12 and 15). This situation is very different from the one in e l g 
perovskites such as LaMn0 3 , where the e g orbitals pda-bond -and thus cause strong 
JT distortion- to different oxygen orbitals than those bonding to A d. Through the 
series of t\ g perovskites, the increasing misalignment of the xy, yz,&nd zx orbitals and 
the theft of O p character by the A ions lead to a decrease of the t 2g bandwidth, W, by 
about 50 % (tables 6 and 8 and figures 10 and 11), and the increasing A d 3z 2^_ 1 and d xy 
character leads to increasing inter-orbital coupling, splitting of the t 2g degeneracy (figure 
14 and equations (17)-(20)), and to the development of a lowest subband (figures 25, 
27 and 29) with a width, Wj, which decreases from W to \W (table 19). To the LDA 
low-energy Hamiltonian defined by the t 2g Wannier functions we have finally added the 
on-site Coulomb interaction terms. 

In order to calculate the correlated spectral densities and, hence, to study the 
Mott transition, Pavarini et al. [21] solved this Hubbard Hamiltonian in the single- 
site dynamical mean-field approximation (DMFT) at temperatures well above those 
where magnetic orderings occur in the titanates. The critical values, U' c , of the average 
Coulomb interaction, U', required to drive the materials into the Mott insulating state 
decrease by a factor two through the series (figure 19 and table 10), thus reflecting not 
only the decreasing W, but also the increasing crystal-field splittings and the decreasing 
width of the lowest subband. For the titanates, the Mott transition occurs essentially 
in the lowest subband, that is, the orbital degeneracy decreases from 3 in SrV0 3 to 
~1 in the titanates (figure 21) and, for U' > U' c , the occupied orbital is the Wannier 
function for the lowest crystal-field level. Using U' ~ 3.7 eV for all four materials, the 
spectral densities in figure 20 not only reproduce the trend that SrV0 3 and CaV0 3 are 
metals with increasing mass enhancements and LaTi0 3 and YTi0 3 Mott insulators with 
increasing gaps, but also the quantitative features of photoemission and BIS spectra. 
In the titanates, the orbital polarization is nearly complete (0.91 in LaTi0 3 and 0.96 
in YTi0 3 ) and the crystal-field splittings and orbital orders agree well with recent 
experiments. The difference between the orbital orders in La and Y titanate (figure 16) is 
quantitative, rather than qualitative, and is caused by the dominating A character being 
^3z? n -i i n LaTi0 3 and being d xy in the more heavily GdFe0 3 -type distorted YTi0 3 . 
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The difference in the onsets of optical conductivities in the titanates is reasonably well 
accounted for (figure 23). The volume reduction needed to make LaTi03 a metal is 
reproduced and that for YTi0 3 is predicted (figure 22). 

In order to calculate the magnetic orderings in the insulators at low temperature 
we have assumed complete orbital order and have applied conventional super-exchange 
theory to our Hubbard Hamiltonians. In accord with experiments, we find LaTi0 3 to 
be a G-type antiferromagnet with fairly isotropic exchange coupling constants, but with 
values three times smaller than those obtained by from spin- wave spectra (table 12). 
We correctly find YTi0 3 to be a ferromagnet with very small coupling constants, but 
their anisotropy is too large (table 11 and figure 24). The reason for the trend from 
antiferro- towards ferromagnetism with increasing GdFeO"3-type distortion is that the 
hopping between the nearest-neighbour Wannier functions for the lowest crystal-field 
level decreases compared with the hoppings to the higher-level Wannier functions. This 
is also the mechanism for creating a lowest subband of reduced width. At the point 
where the super-exchange coupling becomes weak - because it changes sign - the JT 
distortion becomes the controlling factor for the magnetic order, although it hardly 
influences the orbital order. We find that YTi0 3 without JT distortion should be a 
C-type antiferromagnet and that YTi0 3 at 16 GPa pressure, where the JT distortion is 
strongly reduced, should become an A-type antiferromagnet below ~ 100 K (table 16). 

Our results concerning the role of the GdFe0 3 -type and JT distortions for the 
crystal-field splittings, orbital orders, and magnetic couplings are in accord with those 
obtained from studies of model Hamiltonians by Mochizuki and Imada [29]. Both 
theoretical works, together with an increasing number of experimental works, point to a 
crystal-field splitting in the titanates at the order of 200 meV, i.e. much larger that kT 
and the spin-orbit splitting. This puts doubts on the applicability to the titanates of the 
recent orbital-liquid theory [16, 17], which is based on the assumption the t 2g orbitals are 
nearly degenerate. Although the picture we have presented appears consistent from the 
chemical point of view, the fine balance between O-A, O-B, A-B, and B-B covalencies is 
hardly reproduced with sufficient accuracy to describe the magnetic exchange couplings 
well in these t\ g materials. This may be due to the LDA, to our use of the atomic-spheres 
approximation to generate the LDA potentials, and to our NMTO t2 g Wannier functions 
being too extended for use as basis functions for the Hubbard model. Moreover, our 
assumption of complete orbital order for LaTi0 3 may be too crude considering the fact 
that the measured magnetic moment is only 0.57 \ib in this material. Only other types 
of calculations, for instance spin-polarized LDA+DMFT calculations for LaTi0 3 , can 
throw light on this in the future. 

In order to prepare the grounds for future experimental and theoretical clarification, 
and also to ease the applicability of the present theoretical work to the wealth of 
interesting, similar t 2g systems, we have presented a great deal of detail about why 
our results come out the way they do, and we have tabulated the LDA on-site and 
hopping matrix elements (tables 1, 2, 4, and 5). Moreover, in section 3.6 we have given 
an analytical expression in the k + Q representation for the orthorhombic Hamiltonian, 
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into which the self-energy may be included, once it has been extracted from calculations 
or experiments. In the present LDA+DMFT calculations the k-dependence of the 
self-energy was totally neglected, but our success in explaining the Mott transition 
in the series of t\ perovskites indicates that, for these systems, the single-site DMFT 
is indeed a good approximation. On the other hand, since the GdFeC>3-type distortion 
severely influences not only the on-site matrix elements, but also the hopping integrals 
in t 2g systems, this might induce a k-dependence of the self-energy not found in e g 
systems. In the near future it would be useful to get hold of the self-energy matrix for 
real frequencies, because this would give us the correlated bandstructure (k-dependent 
spectral function), which could then be compared with for instance angle- resolved 
photoemission and dHvA experiments. For metals, theoretical work on the electron 
liquid [71] indicates that the self-energy does have a significant k-dependence, and that 
could be checked. In the future it may be possible to use cluster-DMFT calculations to 
evaluate the k-dependence of the self-energy for real systems. 

In this paper we have also demonstrated the use of the new NMTO downfolding 
technique [30] as a tool for generating truly minimal basis sets, Wannier functions in 
particular. It is of course possible to generate sets of Wannier functions which span also 
the oxygen p and transition-ion e g LDA bands as needed when describing high-energy 
excitations. As another extreme, it is sometimes possible to downfold to even fewer than 
three ti g functions and, in such cases, discover tendencies towards symmetry breaking, 
which may be exploited by the Coulomb correlations (figure 25). We may, in general, 
follow the localization process more closely and accurately by generating the Wannier 
functions self-consistently during the cause of a DMFT calculation [51]. 

We are convinced that the rich physics of the materials studied in this paper, as 
well as the computational techniques used to do so, will remain active fields of research 
for years to come. 
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Appendix A. NMTO basis sets 

The iVth-order muffin-tin orbitals (NMTO) method [30] is more intelligible, flexible, and 
accurate than its predecessor, the linear muffin-tin orbitals (LMTO) method [41], a so- 
called fast band-structure method. In the present paper, we use the NMTO method for 
generation of localized Wannier functions for the Kohn-Sham band-structure problem. 

The method constructs a set of local-orbital basis functions which span the solutions 
of Schrodinger's equation -actually, the (scalar) Dirac equation- for a local potential, 
written as a superposition, J2 r vr (|r — R|) , of spherically symmetric potential wells 
with ranges sr, a so-called overlapping muffin-tin potential. This is done by first solving 
the radial equations (numerically) to find y?^/ (e ra , |r — R|) for all angular momenta, /, 
with non- vanishing phase-shifts, for all potential wells, R, and for a chosen set of energies 
spanning the region of interest, e n = €i,....,€n- These energies are the ones shown on 
the right-hand side of figure 4 for the oxygen 2p bands and on the right-hand side of 
figure 10 for the B3d(t 2g ) bands. 

The partial-wave channels, Rim, are now partitioned into active and passive. The 
active ones are those for which one wants to have orbitals in the basis set; i.e., they are 
the chosen one-electron degrees of freedom. For the red bands in figures 3 and 10, the 
B d (t2 9 ) channels are active, while for the red bands in figure 4, the oxygen p channels 
are active. In all figures, the black bands are calculated with a large NMTO set having 
the Op, Bspd, and Aspd channels active. 

For each active channel, Rim, a so-called kinked partial wave, 4>fam (e«, r) , 



is constructed from all the partial waves, <fm (e n , |r — R|) Y\ m \ y — RJ , inside the 
potential-spheres, and from one solution, tyfam ( e n, r ) , of the wave-equation in the 
interstitial, a so-called screened spherical wave. The construction is such that the kinked 
partial wave is a solution of Schrodinger's equation at energy e„ in all space, except at so- 
called hard screening-spheres -which are concentric with the potential-spheres, but have 
no overlap- where it is allowed to have radial kinks in the active channels. It is now clear 
that if we can form a linear combination of such kinked partial waves with the property 
that all kinks cancel, we have found a solution of Schrodinger's equation with energy e n . 
In fact, this kink-cancellation condition leads to the classical method of Korringa, Kohn 
and Rostoker [72] (KKR), but in a general, so-called screened representation and valid 
for overlapping MT potentials to leading order in the potential overlap. The screened 
KKR equations are a set of energy- dependent, homogeneous linear equations, with a 
matrix, Kg^- (e) , whose rows and columns are labelled by the active channels. We 
do not solve this set of secular equations, but proceed a bit differently: 

The major computational task in the screened KKR and NMTO methods is to 
construct the set of envelope functions for the kinked partial waves, the set of screened 
spherical waves, as superpositions of spherical Hankel functions. This is known as the 
real-space calculation of the screened structure-matrix, which is the non-diagonal part 
of the KKR matrix, K (e n ). Now, ipRjm (e n ,r) must join smoothly onto all the passive 
partial waves, z.e.it must have the proper phase shifts. For all the active channels, 
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except the eigenchannel, Rim, it can be forced to vanish at the screening spheres, and 
it is this confinement which makes it localized, provided that this is possible with the 
actual choices of energy, partition between active and passive channels, and screening- 
radii, or. Since the screened spherical wave is required to vanish merely for the other 
active channels, but not for the eigenchannel, it is an impurity solution for the hard- 
sphere solid. In order to obtain maximal localization, the hard spheres are usually 
chosen to be nearly touching. The passive channels are said to be downfolded. 

As an example, let us consider the set of 2p kinked partial waves at the energy e± 
in figure 4 and assume that this set equals the set of NMTOs shown in figures 5 and 6, 
which is approximately true. This set consists of the three p x , p y , and p z orbitals on all 
oxygens in the solid. The 01 p x orbital shown in the 3rd column and upper row of figure 
5 therefore satisfies the following conditions: at its own site, the p y - and p z -projections 
vanish, at all other oxygen sites in the solid, all three oxygen p-projections vanish, 
and all remaining projections, i.e. Osd.., La spdf.., and Ti spd.., are smooth solutions 
of Schrodinger's equation at energy e\. Also the oxygen p-projections are solutions of 
Schrodinger's equation at energy ei, but they have kinks at the screening-spheres. 

Another example are the B d xy , d yZ) and d xz NMTOs shown figure 8. In order to 
generate the orbitals with the convention used in the second row, the active channels 
were specified simply as d xy , d yz , and d xz on each of the four B-sites, with x, y, and z 
referring to the global axes. This is possible when the structure is nearly cubic, because 
then the orbitals can adjust their orientation due to the freedom of, say, the d^-NMTO 
to contain any on-site character, except d yz and d xz , such as for instance d x 2_ y 2 and 
d 3z 2_ 1 . 

The set of NMTOs is formed as a superposition of the kinked-partial-wave sets for 
the energies, e±, 

N 

n=0 Rim 

Note that the size of this NMTO basis set is given by the number of active channels 
and is independent of the number, N + 1, of energy points. The coefficient matrices, 
in equation (A.l) are determined by the condition that the set of NMTOs span 
the solutions, (e i: r) , of Schrodinger's equation with an error 

*\ N) (r) - Hfi ( £i , r) = C W (e« - e ) (e* - Cl ) ... (e* - e N ) (A.2) 

+ o ((ei - e ) (ei - ei) ... (ei - e N )) . 

This is polynomial approximation for the Hilbert space of Schrodinger solutions and 
are the coefficients in the corresponding Lagrange interpolation formula. An NMTO 
with > 0, has no kinks, but merely discontinuities in the (2A r + l)st radial derivatives 
at the screening-spheres for the active channels. 

The Lagrange coefficients, L^, as well as the Hamiltonian and overlap matrices in 
the NMTO basis are expressed solely in terms of the KKR resolvent, K , and its 
first energy derivative, K (e) -1 , evaluated at the energy mesh, e — e±, e^. Variational 
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estimates of the one-electron energies, such as the red bands in figures 3, 4, and 10 
may be obtained from the generalized eigenvalue problem, 

((x iN) \n\x (N) )-^(x (N) \x (N) ))^ = o, 

with 

W = -A + 5>*(|r-R|), 

R 

or as the eigenvalues of the one-electron Hamiltonian matrix, 

H LDA = (x w± \n\x {N)± ) (A.3) 
in the basis of symmetrically orthonormalized NMTOs: 

W N)± ) = \x {N) ) (x w \x {N) y h - (A.4) 

The prefactor, c^ N \ of the leading error of an NMTO set (A. 2) depends on the size 
of the set, and the larger the set, the smaller the prefactor. For that reason we have not 
bothered to indicate the energy mesh used for the black bands in figures 3, 4, and 10. 
Assuming that they are exact, the red bands will touch the black bands quadratically 
at the energy points used to generate the red-band set. That they touch rather than 
cross is by virtue of the variational principle. 

For an isolated set of bands, like the oxygen j9-bands or the B t 2s -bands, and for 
energy meshes spanning the range of those bands, as the number of energy points 
increases and the distances between them decrease, the set of truly minimal NMTOs 
will converge and be exact. Upon orthonormalization, they will therefore form a set of 
localized Wannier functions. 

The construction of a minimal NMTO basis set is different from standard Lowdin 
downfolding. The latter partitions a given, large (say orthonormal) basis into active (A) 
and passive (P) subsets, then finds the downfolded Hamiltonian matrix as: 

(A (e) \H\ A (e)) = (A \H\ A) + (A \H\ P) (P \e - H\ P)' 1 (P \H\ A) , (A.5) 

and finally removes the e-dependence of the downfolded basis by linearizing 
(P \H — e\ P) -1 and treating the term linear in e as an overlap matrix. Obviously, since 
the NMTO set is exact at N + 1 energy points, it is more accurate. Nevertheless, 
since truly minimal NMTOs can be fairly complicated functions, we often use 
equation (A.5) to interpret their (orthonormalized) Hamiltonian in terms of the 
Hamiltonian represented in a larger basis set whose orbitals are simpler and more 
localized. 

For crystals, all calculations except the generation of the screened structure matrix 
are performed in the Bloch k-representation 

x!£fc'> = :^ < A - 6) 

where T labels the L (— > oo) lattice translations and R the active sites in the primitive 
cell. In order to obtain the orbitals and the Hamiltonian in configuration space, Fourier- 
transformation over the Brillouin zone is performed. 
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In the present paper, the orbitals shown are NMTOs before orthonormalization 
(A. 4) because they are (slightly) more localized than the orthonormalized ones. The 
hopping integrals and on-site elements given in the tables are of course matrix elements 
of the orthonormalized Hamiltonian (A. 3). 



Appendix B. Technical details of the LMTO-ASA potential calculations 

Since our present NMTO code is not self-consistent, we used the current Stuttgart TB- 
LMTO-ASA code [73] to generate the LDA potentials. Such a potential in the atomic- 
spheres approximation (ASA) is an overlapping MT-potential, like the one handled by 
the NMTO method, but with the relative overlaps, 

URR ' = iR-R'l ' ( ^ 

limited to about 20%. This limitation comes from the LMTO-ASA+cc method[41], 
which solves Schrodinger's equation by treating the overlap as a perturbation (the 
socalled combined-correction term, cc) and uses screened spherical waves of zero kinetic 
energy in the s-interstitial. Poisson's equation is solved for the output charge density, 
spherically symmetrized inside the same atomic s-spheres. For a given potential, the 
hopping integrals obtained with the NMTO method are more accurate than those 
obtained with the LMTO method, first of all because the NMTOs do not use the zero- 
kinetic-energy approximation in the interstitial region and, secondly, because we use 
N > 1 with well-chosen energy meshes. However, since the on-site matrix elements of 
the LDA Hamiltonian (A. 3) turned out to be crucial for the present study, it is possible 
that our ASA treatment of the potential leads to an -underestimation of the electrostatic 
contribution to the crystal field. 

We now specify our computational set-up. The radii of the potential spheres, sr, 
were dictated by our use of the LMTO-ASA method to generate the LDA potentials. 
In order to limit the overlaps defined by equation (B.l), interstital -or empty- spheres 
(E) were inserted in the non-cubic structures. Table Bl gives the radii of the potential 
spheres. As a result, the overlap between atomic spheres was <16%, between atomic 
and empty spheres <18%, and between empty spheres <20%. We used the guidance 
given by the current version of the code in choosing the potential spheres appropriately. 

With these reasonably large oxygen spheres, the oxygen 2s electrons could be 
treated as part of the core. The self-consistent valence-electron densities were calculated 
with the LMTO bases listed in Table B2. In order to describe properly the A-O-B 
covalency we found it important to downfold the oxygen d partial waves, rather than 
to neglect them (i.e. to approximate them by spherical Bessel functions when solving 
Schrodingers equation, and to neglect them in the charge density). Since the LMTO 
calculations were used to produced the self-consistent charge densities, the energies, em, 
for the linear <pm, 4>ri expansions were chosen at the centres of gravity of the occupied 
parts of the respective DOS i?/-projections. 
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Table Bl. Radii sr of potential spheres in Bohr atomic units. 



AB0 3 


A 


B 


01 


02 


E 


El 


E2 


E3 


E4 


E5 


SrV0 3 


3.97 


2.29 


1.77 
















CaV0 3 


3.34 


2.33 


1.80 


1.80 


1.28 


1.25 


1.17 


1.22 


1.24 




LaTi0 3 [44] 


3.37 


2.51 


1.90 


1.90 


1.50 


1.41 


1.31 


1.29 


1.19 


1.13 


LaTi0 3 [12] 


3.31 


2.52 


1.93 


1.93 


1.61 


1.49 


1.42 


1.28 


1.08 


1.07 


YTi0 3 


2.95 


2.51 


1.92 


1.90 


1.93 


1.83 


1.48 


1.44 







Table B2. LMTO basis sets used in the self-consistent calculation of LDA potential. 
(I) means that the /-partial waves were downfoldcd within in the LMTO-ASA+cc. 





A 


B 01 


02 


E 


El 


E2 


E3 


E4 


E5 


SrV0 3 


s(p)d(f) 


spd (s)p(d) 
















CaV0 3 


s(p)d 


spd (s)p(d) 


(s)p(d) 


s(p) 


s{p) 


s(p) 


s(p) 


s(p) 




LaTi0 3 


s(p)df 


spd (s)p(d) 


(s)p(d) 


s(p) 


s{p) 


s(p) 


s(p) 


s(p) 


s(p) 


YTi0 3 


s(p)d(f) 


spd (s)p(d) 


(s)p(d) 


s(pd) 


s(pd) 


s(p) 


s(p) 







With the benefit of hindsight, we should have inserted the 12 E spheres above the 
octahedron edges also in cubic SrV0 3 . That would have reduced s Sr to 3.46 a.u., a 
value closer to those for the other systems. As a result, the ti g bandwidth for cubic 
SrV0 3 would be reduced by 8%. This seems to be the largest computational "error" of 
the present calculations. This inaccuracy only concerns cubic SrV0 3 and it reduces the 
decrease of rms bandwidth (tables 6 and 8) when going from SrV0 3 to CaV0 3 , from 19 to 
11%. Hence, the trend that CaV0 3 is a more correlated metal than SrV0 3 is somewhat 
diminished. All remaining results, such as those concerning the Mott transition and the 
properties of the titanates, of course remain valid. 

After completion of the calculations, we also found that the optimal trade off 
between the errors caused by the confinement of the O 2s electrons to the sphere and 
by the overlap of spheres, is obtained with a larger oxygen radius, 2.04 a.u.. This leads 
to a small downwards shift of the O 2p band and a 4% decrease of the ti g bandwidth, 
W, for all four materials. Since all materials are influenced the same way, and since 
our value of U was fitted relatively to W, this error has no effect, except that it might 
influence the sensitive values of the exchange-coupling constants calculated in section 6. 

Finally, in the NMTO calculations, the hard-sphere radii, or, for the active channels 
were chosen as 0.7sr. 
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